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PROGRAM  GOALS 


This  grant  is  for  a  twenty-seven  month  investigation  of  the  influence  of 
impurities  and  vortices  on  the  nearby  electronic  structure  in  a  short  coher¬ 
ence  length  superconductor.  The  primary  goal  is  to  develop  self-consistent 
calculations  of  the  electronic  structure  near  various  defects,  focusing  on  the 
local  density  of  states  (LDOS)  measurable  by  an  STM.  These  calculations 
are  to  be  compared  with  analytic  results  in  various  limits.  Pairing  functions 
of  d-wave  and  s-wave  symmetry  are  to  be  considered.  Results  for  the  lo¬ 
cal  electronic  structure  are  expected  to  be  of  great  utility  in  understanding 
macroscopic  measurements  such  as  planar  tunneling  spectra,  critical  currents 
and  the  magnetic  penetration  depth,  but  the  application  of  theoretical  results 
to  these  problems  is  currently  not  funded. 

PROGRESS 

Over  the  past  year  we  have  developed  a  new  Green’s  function  technique, 
which  we  call  the  Koster-Slater  technique,  for  calculating  the  local  electronic 
structure  in  superconductors.  This  technique  allows  us  to  calculate  the 
changes  in  the  order  parameter  self  consistently  and  is  applicable  to  weak 
or  strong-coupling  superconductors  with  arbitrary  order  parameters.  In  a 
practical  numerical  calculation  the  Koster-Slater  technique  allows  for  essen¬ 
tially  arbitrarily  fine  energy  resolution  in  the  calculated  tunneling  spectra, 
and  requires  approximately  one  hour  to  run  for  a  high-temperature  supercon¬ 
ductor  on  a  high-end  desktop  workstation.  In  contrast,  the  Bogoliubov-de 
Gennes  equations  have  lousy  energy  resolution  for  realistic  superconducting 
energy  gaps,  and  we  estimate  would  require  six  months  to  run  for  a  single 
spectrum  (if  there  were  a  computer  with  one  million  MB  of  memory). 

The  initial  application  of  this  technique  was  to  an  isotropically  gapped 
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low-temperature  superconductor  because  of  recent  STM  results  on  the  local 
density  of  states  near  magnetic  impurities  on  niobium  (Don  Bigler’s  group, 
IBM  Almaden).  We  explained  the  qualitative  features  of  the  data  with  a 
calculation  based  on  a  free-electron  model  (attached  reprint).  We  have  also 
explored  many  general  features  of  the  STM  spectra  around  impurities  on 
isotropically  gapped  superconductors  using  the  remarkable  energy  resolution 
of  the  Koster-Slater  technique  (attached  preprint).  These  include  the  spec¬ 
tra  around  non-magnetic  and  magnetic  impurities  in  isotropically  gapped 
superconductors  with  long  and  short  coherence  lengths. 

A  major  theoretical  conclusion  of  this  work  was  that  use  of  the  correct 
detailed  band  structure  of  a  material  is  critical  for  an  accurate  calculation 
of  the  STM  spectra  near  an  impurity.  This  is  the  case  because  the  energy 
scale  of  impurities  is  atomic  (eV),  even  though  the  energy  of  states  created 
within  the  gap  is  of  order  the  gap  (meV).  Hence  the  impurity  potential 
mixes  in  states  of  the  superconducting  host  which  are  very  far  away  from 
the  Fermi  surface,  and  so  a  correct  band-structure  for  those  far  away  states 
is  necessary.  We  are  currently  setting  up  a  realistic  band-structure  model 
for  niobium,  based  on  thirty-four  tight-binding  parameters  going  out  to  the 
third  nearest  neighbors,  which  will  allow  us  to  compare  quantitatively  with 
the  experimental  data. 

We  have  also  produced  an  analytic  model  for  the  impurity  within  the 
s-wave  superconductor  based  on  a  delta-function  potential.  Previous  models 
assumed  particle-hole  symmetry,  and  thus  produced  pathological  results  for 
the  local  STM  spectra.  Our  new  analytic  model  has  a  single  parameter 
which  accounts  for  the  actual  case  of  no  particle-hole  symmetry.  Thus  the 
qualitative  behavior  of  the  local  STM  spectra  can  be  understood  without 
a  detailed  numerical  calculation.  This  model  is  detailed  in  the  attached 
preprint. 

V 

Once  the  utility  of  the  new  technique  was  demonstrated  on  isotropically 
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gapped  superconductors  we  applied  it  to  anisotropically  gapped  supercon¬ 
ductors  with  a  tight-binding  band  structure.  The  order  parameters  chosen 
were  anisotropic  s-wave  and  d-wave.  STM  spectra  for  magnetic  and  non¬ 
magnetic  impurities  in  these  superconductors  have  been  calculated,  and  pa¬ 
pers  on  these  results  are  in  preparation.  For  magnetic  impurities  we  find  that 
the  sharp  localized  states  within  the  gap  which  are  present  for  isotropically 
gapped  superconductors  broaden  into  resonant  states  in  the  anisotropically 
gapped  superconductors.  These  resonant  states  have  highly  asymmetric  line- 
shapes  and  are  present  for  both  d-wave  and  anisotropic  s-wave  order  param¬ 
eters.  Resonant  states  also  are  present  for  non-magnetic  impurities  for  both 
types  of  order  parameter. 

The  resonant  states  that  we  find  should  influence  planar  tunneling  spectra 
when  the  impurity  concentration  is  high  enough.  Since  YBCO  typically  has 
a  disordered  surface,  it  is  possible  that  these  resonant  states  may  be  a  cause 
of  the  two-gap  features  seen  in  many  tunneling  experiments.  This  type  of 
structure  in  the  density  of  states  of  the  dirty  superconductor  should  produce 
signatures  in  an  assortment  of  macroscopic  quantities,  such  as  the  specific 
heat  and  the  magnetic  penetration  depth.  In  general,  features  which  appear 
as  double  gaps  may  originate  from  an  anisotropic  order  parameter  (s  or  d) 
in  the  presence  of  impurities. 

After  submitting  our  results  on  anisotropically  gapped  superconductors 
we  plan  to  investigate  the  other  types  of  defects  discussed  in  our  proposal: 
step  edges  and  vortices.  An  investigation  of  the  step  edge  will  be  straight¬ 
forward  now  that  the  single  impurity  problem  is  under  control.  Treating  the 
twisting  of  the  phase  of  the  order  parameter  around  a  vortex  core,  however, 
will  require  a  modification  of  the  Koster-Slater  technique.  We  expect  this 
problem  to  be  solvable  without  great  difficulty. 


4 


LIST  OF  ATTACHED  PAPERS 


1.  Michael  E.  Flatte  and  Jeff  M.  Byers,  “Local  Electronic  Structure  of  a 
Single  Magnetic  Impurity  in  a  Superconductor” ,  Physical  Review  Let¬ 
ters  78,  3761-3764  (1997).  Development  of  the  Koster-Slater  technique 
and  first  application  (magnetic  impurities  on  niobium). 

2.  Michael  E.  Flatte  and  Jeff  M.  Byers,  “Local  Electronic  Structure  of  De¬ 
fects  in  Superconductors” ,  Physical  Review  B,  in  press.  Further  devel¬ 
opment  of  the  Koster-Slater  technique  and  application  to  magnetic  and 
non-magnetic  impurities  in  isotropically  gapped  superconductors  with 
long  and  short  coherence  lengths.  Outline  of  technique  for  anisotropi- 
cally  gapped  superconductors. 

SEMINARS  AND  CONFERENCE  TALKS 

1.  “Local  Properties  of  Impurities  in  Superconductors”,  Solid  State  Sem¬ 
inar  at  Iowa  State  University,  October  3,  1996. 

2.  “Local  Properties  of  Impurities  in  Superconductors” ,  Michael  E.  Flatte 
and  Jeff  M.  Byers,  American  Physical  Society  March  Meeting,  Kansas 
City,  March  21,  1997. 

3.  “Local  Properties  of  Impurities  in  Superconductors” ,  Solid  State  Sem¬ 
inar  at  Harvard  University,  April  24,  1997. 

INVITED  PAPERS 

1.  I  have  been  invited  to  write  a  review  article  entitled  “Local  Proper¬ 
ties  of  Defects  in  Superconductors”  for  the  Solid  State  Physics  series, 
published  by  Academic  Press. 


5 


Volume  78,  Number  19 


P  H  Y  S I C  A  C  R  E  V I E  W  LETTER  S 


12  May  1997 


Local  Electronic  Structure  of  a  Single  Magnetic  Impurity  in  a  Superconductor 

Michael  E.  Flatte*  and  Jeff  M.  Byers^ 
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^Materials  Physics,  Naval  Research  Laboratory,  Washington,  D.C.  20375 
(Received  11  July  1996) 

We  present  the  first  three-dimensional,  self-consistent  calculation  of  the  electronic  structure  near 
a  strong  classical  magnetic  impurity  in  a  superconductor.  Localized  excited  states  are  found  within 
the  energy  gap  which  are  half  electron  and  half  hole.  The  spatial  structure  of  the  positive- frequency 
(electronlike)  spectral  weight  (or  local  density  of  states)  can  differ  strongly  from  that  of  the  negative 
frequency  (holelike)  spectral  weight.  The  effect  of  the  impurity  on  the  continuum  states  above  the 
energy  gap  is  calculated  with  good  spectral  resolution  for  the  first  time.  For  sufficiently  strong  impurity 
potentials,  the  order  parameter  may  change  sign  at  the  impurity  site.  [80031-9007(97)03084-6] 

PACS  numbers:  74.80.-g,  74.25. Ha,  61.16.Ch 


Magnetic  impurities  have  a  dramatic  effect  on  super¬ 
conductivity.  Most  works  on  the  experimental  and  theo¬ 
retical  effects  of  magnetic  impurities  have  focused  on  bulk 
thermodynamic  quantities,  such  as  the  reduction  of  the 
superconducting  transition  temperature  Tc  with  increas¬ 
ing  magnetic  impurity  concentration  [1].  Theoretical  ap¬ 
proaches  to  treat  the  thermodynamic  effects  of  impurities 
include  Born  scattering  (Abrikosov-Gor’kov  theory)  [2] 
and  approximate  solutions  to  all  orders  in  the  impurity  po¬ 
tential  for  classical  [3-5]  and  quantum  [6,7]  spins.  A  key 
issue  addressed  in  Ref.  [3]  was  the  evolution  of  localized 
excited  states  (not  obtainable  in  Bom  theory)  around  impu¬ 
rities  into  an  impurity  band,  eventually  leading  to  gapless 
superconductivity  at  high  enough  impurity  concentrations. 
The  effects  of  impurities  on  bulk  properties  have  also  been 
treated  within  a  strong-coupling  formalism  (e.g..  Ref.  [8]), 
but  only  very  recently  self-consistently  and  beyond  the 
Bom  approximation  [7].  Concern  about  bulk  properties 
in  the  above  and  related  work  did  not  extend  to  properties 
very  near  the  impurity. 

Among  the  first  properties  calculated  in  the  vicinity  of  an 
impurity  were  the  structures  of  screening  clouds  around  a 
charged  impurity  [9,10]  and  a  magnetic  impurity  [10,1 1]  in 
a  superconductor  (characterized  by  exponentially  decaying 
Friedel-like  oscillations).  The  oscillation  of  the  order 
parameter  around  a  magnetic  impurity  was  first  evalu¬ 
ated  without  self-consistency  [12,13].  A  self-consistent 
calculation  of  the  order  parameter  at  the  impurity  and 
very  far  away  for  weak  impurity  potentials  was  done  by 
Schlottmann  [14].  Interest  in  local  properties  near  impu¬ 
rities  has  been  revived  by  advances  in  scanning  tunneling 
microscopy  (STM)  near  impurities  embedded  in  a  metallic 


medium  [15]. 

The  differential  conductivity  measured  by  a  STM  at  a 
point  X,  voltage  V,  and  temperature  T  can  be  related  to 
the  local  density  of  states  (LDOS)  at  the  tip  location  as 
follows: 


dI(x,V,T) 

dV 


dco  dn{co) 
IttNo  3(0 


Im  Gcr(x,x;  a)). 


(1) 


Here  a  is  the  spin  direction  of  the  electron  (±1/2), 
No  is  the  normal-state  density  of  states  for  each  spin 
direction  at  the  Fermi  level,  and  dljdY  is  in  units  of  the 
normal-metal’s  differential  conductivity.  n{o})  =  [1  + 
exp -  eV)y^  where  e  is  the  electron’s  charge  and 
P  =  1/r.  The  LDOS  is  proportional  to  the  imaginary 
part  of  the  retarded  Green’s  function  fully  dressed  by  the 
interaction  of  the  electronic  system  with  the  impurity.  The 
LDOS  was  calculated  near  impurities  for  superconductors 
with  isotropic  and  anisotropic  gaps  within  the  Bom  ap¬ 
proximation  [16].  A  self-consistent  calculation  of  the 
LDOS  of  a  vortex  core  in  a  superconductor  was  performed 
using  the  Schrodinger-like  Bogoliubov-de  Gennes  (BdG) 
equations  [17,18]  and  compared  with  measurements  by  an 
STM  on  superconducting  NbSe2  [19]. 

Recent  STM  results  near  a  magnetic  impurity  [20]  indi¬ 
cate  that  (1)  discrete  states  are  evident  within  the  energy 
gap  and  the  LDOS  associated  with  them  is  asymmetric 
with  voltage  around  V  =  0,  and  (2)  the  LDOS  becomes 
indistinguishable  from  the  bulk  density  of  states  within  a 
distance  greater  than  the  Fermi  wavelength,  but  much  less 
than  the  coherence  length. 

We  will  present  a  fully  self-consistent  calculation  of 
the  local  electronic  structure  near  a  magnetic  impurity, 
which  is  based  on  a  Koster-Slater-like  Green’s  function 
technique  rather  than  the  BdG  equations  [21].  Since 
the  original  application  of  the  Green’s  function  algorithm 
to  localized  vibrational  modes  [22],  it  has  been  applied 
to  numerous  problems,  including  conduction  electrons  in 
metals  [23]  and  impurity  states  in  magnets  [24],  but  not 
to  superconductors.  The  Green’s  function  method  works 
within  a  sphere  whose  radius  is  determined  by  the  range 
of  the  inhomogeneous  potential.  In  contrast,  the  BdG 
equations  need  to  be  solved  in  a  large  sphere  whose  radius 
is  determined  by  the  required  spectral  resolution.  In  es¬ 
sence  we  shall  invert  in  real  space  the  Dyson  equation  [25], 

J  dx"[S{\  -  x")  -  go.(x,x",6>)V(x")]G<^(x",x';w) 
=  go.(x,  x';  (o) ,  (2) 
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where 


go.(x,x';a> 


g^(x,x';  w) 
/;(x,x';  -(o) 


(3) 


the  same  relationship  exists  among  Go-,  Ga,  and  F^,  and 


V(x)  =  ( 


o-V5(x)  +  Vo(x)  5A(x) 

5A(x)  crV5(x)  -  Vo(x) 


(4) 


5A(x)  =  A(x)  -  Ao,  where  A(x)  is  the  position- 
dependent  order  parameter  and  Ao  is  the  order  parameter 
of  the  homogeneous  superconductor.  VjCx)  is  a  localized, 
spin-dependent  potential,  such  as  one  originating  from  a 
classical  spin.  Vo(x)  is  a  localized  nonmagnetic  potential. 
The  quantization  direction  of  the  superconductor’s  spins 
is  chosen  parallel  to  the  direction  of  the  classical  spin 
and  so  that  spin-up  electrons  are  attracted  to  the  impurity. 

+  G -a  due  to  the  differences  in  the  potential  in 
Eq.  (4).  For  a  BCS  superconductor  with  an  isotropic 
order  parameter  in  a  parabolic  band,  for  <y  much  smaller 
than  the  Fermi  energy,  the  homogeneous  retarded  Green’s 
functions  are 


go-(x,x';<M)  =  - 


'•  =  - 


X  ( cos  kfr  + 
TTNot^o 


0) 


sin  kpr 


)■ 


X  sin  kpr , 


(5) 


where  r  =  |x  -  x'|,  Ao  is  the  order  parameter,  and  ^  = 
hvpliTh.0  is  the  coherence  length.  These  expressions  are 
valid  for  a  above  and  below  Ao  so  long  as  the  imaginary 
parts  of  both  fa-  and  are  multiplied  by  sgn  <«;.  For  a 
pairing  strength  y,  the  self-consistency  equation  is 

/oo 

d(UM(w)ImFo-(x,x;<w)  —  Ao  .  (6) 

-oo 

One  strength  of  this  formalism  is  its  reliance  on 
the  short-range  nature  of  the  inhomogeneous  potential. 
The  solution  of  Eq.  (2)  requires  inverting  the  frequency- 
dependent  real-space  matrices  Mo-(£u)  =  5(x  —  x')  — 
go-(x,  x';  (y)V(x')-  For  x'’s  where  V(xO  is  negligible,  that 
portion  of  Mo-(<w)  is  triangular  and  trivially  invertible. 
Hence  the  numerical  complexity  of  inverting  Mo-(a))  is 
determined  by  the  range  of  V(x),  which  is  governed  by 
the  range  of  the  longest-range  potential,  which  in  this 
Letter  will  be  SA(x). 

As  a  special  case  of  Eq.  (2),  the  energies  of  the  localized 
states  within  the  gap  correspond  to  those  aj  >  0  where 
detMo-(<u)  =  0.  We  model  the  impurity  with  a  Gaussian 
Vs(x)  of  range  Figure  1  shows  the  dependence  of  the 
localized  state  energies  for  the  first  two  angular  momentum 
channels  on  the  strength  of  the  magnetic  potential,  — 
-ttNoIct  f  dxVs(x)l  We  note  that  a  localized  state  implies 
local  pair-breaking  excitations  of  energy  less  than  2Ao.  If 


FIG.  1.  Localized  excited  state  poles  for  angular  momentum 
channels  €  =  0  and  1  (solid  lines).  Dashed  lines  show  results 
of  Ref.  [3].  At  Vs  =  Vso,  the  character  of  the  €  —  0  state 
changes  from  spin  up  to  spin  down. 


detMo-(«)  =  0,  then  detM-o-C-w)  =  0,  so  the  poles  of 
Ga-icv)  must  come  in  symmetric  pairs  (cr  =  ±  1/2)  around 
(0  =  0.  The  quasiparticle  state  for  small  Vs  corresponding 
to  these  poles  consists  of  an  electron  (a»  >  0)  in  the  spin- 
up  (t)  band  (the  band  attracted  to  the  classical  spin)  and 
a  hole  (co  <  0)  in  the  spin-down  (i)  band.  A  hole  in  the 
down  band  has  spin  up  and  is  attracted  to  the  classical  spin. 
Spin  is  a  good  quantum  number  for  the  quasiparticle,  which 
is  spin  up.  The  spatially  integrated  spectral  weight  is  one- 
half  each  for  the  electron  and  hole  components. 

As  the  potential  strength  increases,  the  excitation  en¬ 
ergy  of  each  angular  momentum  state  €  decreases,  and 
then  at  some  critical  value  {v*()  begins  to  increase.  For 
Vs  >  v*  the  ground  state  contains  a  spin-up  quasiparticle 
bound  to  the  classical  spin  [5].  The  low-energy  excitation 
now  corresponds  to  removing  that  quasiparticle,  which 
is  a  spin-down  excitation.  This  qualitative  behavior  of 
the  excitation  energies  can  be  extracted  from  the  Shiba 
model  [3],  where  the  magnetic  potential  is  modeled  by  a 
delta  function  at  the  impurity  site.  In  Ref.  [3]  the  zeros 
of  Mo-M  are  found,  neglecting  the  frequency  symmet¬ 
ric  component  (whose  significance  is  discussed  below)  of 
Re gix,  x';  <u).  These  poles  are  also  shown  in  Fig.  1. 

The  decrease  in  the  quasiparticle  energy  by  the  classical 
spin  is  proportional  to  the  overlap  of  the  atomic-scale 
potential  and  the  impurity  state.  The  state  is,  aside  from  an 
r~^  falloff,  extended  over  a  volume  whose  linear  scale  is 
the  coherence  length.  Thus  the  overlap,  and  consequently 
the  quasiparticle’s  energy  reduction,  is  proportional  to 
VsAo-  When  this  energy  exceeds  the  energy  required  to 
create  a  quasiparticle  (A<,),  the  ground  state  contains  a 
quasiparticle.  Hence  v*  is  approximately  independent  of 
Ao  [3],  and  is  determined  primarily  by  the  normal-state 
properties  of  the  superconductor. 

In  the  normal  state  the  spectral  weight  of  the  spin- 
up  band  at  the  Fermi  energy  has  a  peak  at  the  origin, 
while  that  of  the  spin-down  band  is  pushed  away  from  the 
origin.  In  the  superconductor  this  asymmetry  is  evident 
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in  the  quasiparticle  spectral  weights.  Figure  2  shows  the 
calculated  angular  momentum  €  =  0  spectral  weights  in 
the  up  and  down  bands  for  two  values  of  Vs,  indicating 
that  the  asymmetry  of  the  spatial  structure  of  the  spectral 
weights  becomes  more  pronounced  with  increasing  v^. 
Also  shown  are  the  €  =  1  spectral  weights  for  = 
0.875.  The  asymmetric  localized-state  spectral  weights 
produce  differential  conductivities  which  are  asymmetric 
as  well.  Figure  3  shows  differential  conductivities  for 
Vs  =  7/8  as  a  function  ofthe  voltage  and  the  distance  from 
the  impurity.  By  r  =  5kp^  the  spectrum  has  recovered 
to  the  homogeneous  (bulk)  spectrum.  The  differential 
conductivity  farther  away  can  be  easily  recovered  by 
constructing  the  self-consistent  T  matrix,  Y{x')S{\'  — 
x")  +  V(x')G(x',x";  (i))V(x"),  directly  from  the  dressed 
Green’s  functions. 

The  asymmetric  behavior  of  the  electronlike  and  hole¬ 
like  spectral  weights  due  to  the  difference  between  the 
spin-up  and  spin-down  band  states  does  not  emerge  from 
the  model  of  Ref.  [3].  The  Shiba  model,  due  to  the  as¬ 
sumption  of  particle-hole  symmetry  implied  by  setting  the 
symmetric  real  part  of  g(0, 0;  <w)  to  zero,  predicts  identical 
spatial  dependences  of  the  spectral  weights  for  the  up  and 
down  bands.  As  pointed  out  by  Roster  and  Slater  [23], 
proper  treatment  of  the  symmetric  real  part  of  g(0,0;  w) 
is  essential  for  obtaining  the  local  electronic  properties 
around  an  impurity.  The  quantity  Reg(0,0;  (o)  is  simply 
the  Kramers-Kronig  transform  of  the  density  of  states. 
The  spatially  dependent  form  of  Reg(x,  x';  co)  determines 
how  the  continuum  states  are  distorted  by  the  impurity 
potential  to  form  the  bound  states.  Neglecting  the  sym¬ 
metric  part  of  Re  g(x,  x';  co)  corresponds  to  a  particle-hole 
symmetric  density  of  states  around  the  Fermi  energy.  The 
potentials  considered  here  are  locally  strong,  and  mix  in 
states  of  the  homogeneous  superconductor  with  energies 
much  farther  from  the  Fermi  surface  than  Ao  •  Often  a  be¬ 
nign  approximation  when  bulk  properties  are  being  calcu- 


FIG.  2.  Spectral  weighjts  for  the  €  =  0  state’s  up-band  com¬ 
ponent  (electronlike  for  Vs  <  v*)  and  down-band  component 
(holelike  for  v,  <  v*)  as  a  function  of  position  for  two  values 
of  the  magnetic  potential.  Also  shown  are  the  spectral  weights 
for  the  up  and  down  bands  for  f  =  1  for  —  0.875. 


lated,  particle-hole  symmetry  leads  to  pathological  prop¬ 
erties  when  calculating  local  electronic  structure. 

We  introduce  a  new  parameter  into  the  Shiba  model, 
a,  which  is  the  <u-symmetric  part  of  Re(g(0,0;<u  = 
0))/'7rNo.  The  brackets  indicate  that  g(0,x;tu  =  0)  has 
been  averaged  over  a  small  volume  given  by  the  range  of 
the  potential  modeled  by  the  delta  function.  For  a  A  0 
the  ratio  for  the  localized  quasiparticle  of  the  spectral 
weight  of  the  up  band  at  the  impurity,  At(0),  to  the  spectral 
weight  of  the  down  band  at  the  impurity,  Aj(0),  is 

A|(0)  ^  1  -I-  2avs  +  (1  4-  a^)v^ 

Ai(0)  1  -  2aVs  +  (1  +  a'^)v^ 

The  introduction  of  a  does  not  change  the  localized  state 
energies  qualitatively. 

We  now  discuss  the  structure  of  A(x).  Figure  4  shows 
5A(x)’s  for  two  values  of  Uj.  While  oscillating  with 
wavelength  —Trkf',  (5A(x)  falls  off  to  a  negligible 
potential  within  lOA:^ '.  This  justifies  performing  the 
numerical  calculation  in  a  sphere  of  that  radius.  A  typical 
radial  grid  of  100  points  provides  a  numerically  robust 
solution.  The  self-consistent  solution  also  depends  on 
the  value  of  No^o/k},  which  for  niobium  is  3.6  X  10“^ 
[26].  Relaxation  of  the  order  parameter  on  the  atomic 
scale  has  also  been  seen  in  square-lattice  calculations  of 
the  BdG  equations  for  nonmagnetic  impurities  in  quasi- 
two-dimensional  d-wave  superconductors  [27],  which 
should  behave  in  a  somewhat  similar  fashion  to  magnetic 
impurities  in  5-wave  superconductors. 

As  shown  in  Fig.  4,  for  large  values  of  Vs,  A(x  = 
0)  <  0.  Sign  changes  in  A,  which  could  be  seen  in  pair 
tunneling,  have  been  suggested  for  magnetic  impurities  in 


FIG.  3.  Differential  conductivity  (dl/dV),  relative  to  the 
normal  metal  as  a  function  of  the  voltage  (in  units  of  A„)  and 
the  distance  from  the  impurity  (in  units  of  k/r  )  for  Vj  =  7/8. 
The  presence  of  the  large  peak  near  the  origin  on  the  positive- 
frequency  side  of  the  spectrum  indicates  that^the  impurity  is  not 
stronger  than  the  critical  strength  (Uj  <  Vso)-  P  =  7.5/A„, 
which  for  niobium  corresponds  to  T  ~  2  K. 
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FIG.  4.  Change  in  the  local  order  parameter,  5A(x),  for 
two  values  of  the  magnetic  potential.  The  change  becomes 
negligible  beyond  10A:^‘.  Inset:  Order  parameter  at  the 
impurity,  A(0),  as  a  function  of  indicating  a  discontinuous 
change  at  v^q. 


the  barriers  of  Josephson  junctions  [28].  Our  sign  change 
in  A(0)  occurs  (at  T  =  0)  precisely  at  t;*o-  The  symmetry 
of  Eqs.  (2)-(6)  under  co  — ►  -co  and  ^  implies  that 

ImFt(r,r,6;)  =  -~ImF|(r,  r, -^y).  As  the  pole  in  the 
spin-up  band  goes  from  electronlike  (co  >  0)  to  holelike 
(oj  <  0)  and  the  pole  in  the  spin-down  band  goes  from 
holelike  to  electronlike,  the  contribution  to  A(0)  from  the 
localized  state  changes  sign  abruptly.  A(0)  as  a  function 
of  Vs  is  shown  in  the  inset  of  Fig.  4.  At  F  >  0  the 
transition  would  be  smoothed. 

The  behavior  of  A(0)  as  a  function  of  Vs  comes 
from  the  introduction  at  of  a  quasiparticle  into  the 
ground  state  of  the  system.  The  spin-up  quasiparticle 
localized  near  the  impurity  in  the  ground  state  suppresses 
the  local  order  parameter.  Exciting  the  low-energy  state 
(removing  the  spin-up  quasiparticle)  causes  a  pair  to  enter 
the  condensate.  We  find  that  exciting  the  low-energy 
state  for  Vs  >  VsO  increases  A(0)  and  reduces  the  spin 
of  the  superconductor  near  the  impurity.  Typically  the 
excitation  of  quasiparticles  reduces  A(x)  (which  is  the 
case  for  Vs  <  Vs0\ 

We  have  calculated  the  electronic  properties  near  a 
magnetic  impurity  embedded  in  a  superconductor,  and 
used  these  results  to  explain  STM  measurements  around 
a  magnetic  atom  on  a  surface.  The  asymmetry  of  the 
peak  heights  can  arise  from  particle-hole  asymmetry  in 
the  normal  state  of  the  superconductor.  The  recovery  of 
the  spectrum  to  its  bulk  structure  within  a  few  atomic 
spacings  is  due  to  the  decay  of  the  impurity  state 
spectral  weight,  and  the  change  in  sign  of  the  order 
parameter  for  sufficient  Vs  occurs  because  a  quasiparticle 
is  bound  in  the  ground,  state.  None  of  these  aspects  of 
the  system  change  qualitatively  at  a  surface.  Thus  the 
calculations  presented  here  should  have  great  value  for 
the  experimental  (surface)  configuration. 
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The  electronic  structure  near  defects  (such  as  impurities)  in  superconductors  is  explored  using  a 
new,  fully  self-consistent  technique.  This  technique  exploits  the  short-range  nature  of  the  impurity 
potential  and  the  induced  change  in  the  superconducting  order  parameter  to  calculate  features  in 
the  electronic  structure  down  to  the  atomic  scale  with  unprecedented  spectral  resolution.  Magnetic 
and  non-magnetic  static  impurity  potentials  are  considered,  as  well  as  local  alterations  in  the  pairing 
interaction.  Extensions  to  strong-coupling  superconductors  and  superconductors  with  anisotropic 
order  parameters  are  formulated. 


L  INTRODUCTION 

Low-temperature  superconductors  almost  always  have 
a  high  concentration  of  non-magnetic  impurities.  Even 
in  the  dirty  limit,  however,  where  the  mean  free  path 
is  shorter  than  the  coherence  length,  superconductivity 
enduresL  This  phenomenon  can  be  understood  by  gener¬ 
alizing  BCS  pairing^  to  be  between  degenerate  Kramers 
partnered  states  in  a  time-reversal  invariant  system^. 
Magnetic  impurities,  which  do  break  time-reversal  invari¬ 
ance,  have  more  profound  effects  on  the  superconductor 
in  dilute  concentrations  than  non-magnetic  impurities, 
lowering  the  critical  temperature  Tc'^^®  and  producing 
localized  states  within  the  gap^”^^  which  at  sufficient 
concentrations  hybridize  to  produce  gaplessness^^.  In 
the  course  of  the  investigation  of  the  effects  of  impuri¬ 
ties  on  superconductivity  during  the  last  four  decades, 
the  primary  emphasis  has  been  the  influence  of  impuri¬ 
ties  on  bulk  properties.  These  effects  have  been  treated 
within  a  strong- coupling  formalism  (e.g.  Ref.  15),  but 
only  very  recently  self-consistently  and  beyond  the  Born 
approximation^®.  The  above  work  was  primarily  con¬ 
cerned  with  bulk  or  impurity- averaged  characteristics 
and  ignored  the  spatial  structure  of  electronic  properties 
near  to  the  impurity. 

Among  the  first  local  properties  calculated  in  the  vicin¬ 
ity  of  an  impurity  in  a  superconductor  were  the  structures 
of  screening  clouds  around  a  charged  impurity^^’^®  and 
a  magnetic  impurity^®’^®  in  a  superconductor  (charac¬ 
terized  by  exponentially-decaying  Priedel-like^®  oscilla¬ 
tions).  The  oscillation  of  the  order  parameter  around 
a  magnetic  impurity  was  first  evaluated  without  self- 
consistency^^"*^®.  A  self-consistent  calculation  of  the  or¬ 
der  parameter  at  the  impurity  and  very  far  away  for  weak 
impurity  potentials  was  reported  by  Srhlottmann^"^  for 
magnetic  impurities  and  Ferter^®  for  non-magnetic  im¬ 
purities. 

Interest  in  local  properties  near  impurities  in  super¬ 
conductors  has  been  revived  by  the  capability  of  scan¬ 
ning  tunneling  microscopy  (STM)  to  perform  localized 


spectroscopic  measurements.  The  spatially-dependent 
differential  conductivity  around  a  non-magnetic  impu¬ 
rity  at  the  surface  of  a  superconductor  has  been  con¬ 
sidered  theoretically^®  for  both  isotropic  and  anisotropic 
order  parameters.  This  differential  conductivity  is  pro¬ 
portional  to  the  local  density  of  states  (LDOS)  around 
the  impurity.  In  Ref.  26  the  impurity  was  modeled 
as  a  point  defect,  and  spatial  oscillations  in  the  LDOS 
at  various  voltages  were  calculated.  These  oscillations 
are  the  superconductor’s  analog  of  oscillations  in  the 
LDOS  created  by  step  edges  and  defects  on  noble  metal 
surfaces^^’^®.  The  superconductor’s  LDOS  oscillations 
would  allow  one  to  measure  the  anisotropy  of  the  su¬ 
perconductor’s  order  parameter.  The  conductance  oscil¬ 
lations  for  voltages  just  above  a  gap  minimum  or  maxi¬ 
mum  are  strongly  pronounced  in  the  real-space  directions 
corresponding  to  the  momenta  of  the  gap  minimum  and 
maximum.  Calculations  followed  which  considered  sharp 
energy  features  in  the  scattering  process,  such  as  reso¬ 
nant  states^®’®®.  Nevertheless,  an  important  assumption 
of  these  calculations  has  remained  unchecked  in  detail, 
that  the  electronic  distortions  induced  by  the  impurity 
are  local,  including  the  deformation  of  the  order  parame¬ 
ter.  Self-consistent  calculations  using  the  Bogoliubov-de 
Gennes  (BdG)  equations®^  followed  for  two-dimensional 
systems®^”®®,  but  have  been  hampered  by  finite-size  ef¬ 
fects. 

A  magnetic  impurity  differs  from  a  non-magnetic  im¬ 
purity  in  that  a  localized  state  exists  around  it'^~^®’®®. 
The  first  calculations  of  the  LDOS  of  the  localized  and 
continuum  states  around  a  magnetic  impurity  were  per¬ 
formed  recently  both  with  a  simplified  analytic  model 
and  numerically  via  the  new  self-consistent  technique®"^ 
described  in  this  paper.  Calculations  of  the  LDOS  of  the 
localized  state  with  angular  momentum  quantum  num¬ 
ber  ^  =  0  were  performed  with  a  slightly  different  ana¬ 
lytic  model  and  numerically  in  two  dimensions  via  the 
self-consistent  BdG  equations®®.  These  calculations  were 
motivated  by  preliminary  experimental  results  around  a 
Gd  atom  on  a  niobium  surface®®. 
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Another  type  of  defect  is  a  local  change  in  the  pair¬ 
ing  interaction.  The  resulting  spatially-dependent  order 
parameter  can  then  distort  the  LDOS.  An  order  param¬ 
eter  suppression  can  even  localize  states,  as  in  a  vortex 
core'*^’^^ 

The  local  electronic  properties  of  all  of  these  de¬ 
fects  can  be  calculated  self-consistently  from  the  Gor’kov 
equation^^  without  further  approximation  with  the  new 
technique  introduced  in  Ref.  37.  Our  technique  for  calcu¬ 
lating  the  electronic  structure  around  a  defect  in  a  super¬ 
conductor  is  related  to  the  Koster-Slater  inversion  tech¬ 
niques  for  determining  the  local  electronic  structure  of 
impurities  in  metals^^’^^.  Since  its  original  application 
to  localized  vibrational  modes^^,  this  algorithm  has  been 
applied  to  numerous  problems  including  deep  levels  in 
semiconductors^^  and  impurity  states  in  magnets"^^.  The 
new  Koster-Slater  technique  separates  space  around  the 
defect  into  two  regions:  the  near  field  and  the  far  field. 
The  far  field  is  a  region  distant  enough  from  the  defect 
that  the  potential  is  insignificant  and  the  order  parame¬ 
ter  has  relaxed  back  to  its  homogeneous  value.  The  near 
field  is  the  region  close  to  the  defect  where  the  potential  is 
finite  or  the  order  parameter  is  distorted.  In  essence,  the 
Gor’kov  equation  that  determines  the  Green’s  functions 
of  the  inhomogeneous  superconductor  is  inverted  in  the 
real  space  region  of  the  near  field.  This  paper  describes 
the  technique  in  detail,  expands  on  an  analytic  model 
introduced  in  Ref.  37  and  reports  several  calculations  of 
the  properties  of  the  defects  described  above. 

Previous  attempts  to  calculate  the  structure  of  the  near 
field  have  often  used  the  BdG  equations^^  Other  theories 
of  inhomogeneities,  such  as  Ginzburg- Landau  theory^^  or 
the  Eilenberger  equations^®,  treat  the  spatial  degrees  of 
freedom  as  coarse-grained  over  the  superconductor’s  co¬ 
herence  length.  Coarse-grained  approximations  are  not 
appropriate  for  considering  electronic  structure  on  the 
atomic  scale  near  a  defect.  The  BdG  equations  are  gen¬ 
eralized  Schrodinger  equations  for  the  electron  and  hole 
wavefunctions  of  a  quasiparticle,  and  are  valid  for  a  su¬ 
perconductor  with  an  arbitrarily-varying  order  parame¬ 
ter,  only  constrained  by  the  validity  of  BCS  theory. 

Unfortunately,  these  equations  have  significant  prac¬ 
tical  difficulties  as  well.  Despite  qualitative  success 
modeling  STM  measurements  of  a  single  vortex  in  su¬ 
perconducting  NbSe2^^,  calculations  of  the  electronic 
structure^ using  the  BdG  equations  are  hampered  by 
the  difference  in  energy  scales  between  the  Fermi  energy 
and  the  order  parameter.  Since  the  BdG  equations  are 
solved  numerically  for  a  finite  system,  the  difficulty  of 
the  calculation  is  determined  by  the  necessary  spectral 
resolution.  The  key  energy  scale  which  must  be  resolved 
is  the  superconducting  gap.  Thus,  the  numerical  diffi¬ 
culty  increases  as  the  ratio  pf  the  Fermi  energy  to  the 
gap  becomes  large.  Hence  the  band  structure  assumed 
for  the  superconductor  must  be  somewhat  unrealistic  (for 
Refs.  54,  55  the  Fermi  wavelength  was  approximately 
lOOA,  which  is  inappropriately  large^®).  This  limitation 
extends  to  calculations  of  the  interaction  between  a  vor¬ 


tex  and  an  impurity^^,  the  characteristics  of  the  vortex 
lattice^^’^®,  and  work  on  a  non-magnetic  impurity^^”^^ 
and  a  magnetic  impurity^®  in  a  two-dimensional  s-wave 
or  d-wave  superconductor.  In  contrast,  the  computa¬ 
tional  requirements  of  the  Koster-Slater  technique  are 
determined  by  the  range  of  the  impurity  potential,  rather 
than  the  necessary  spectral  precision. 

In  the  Section  II  of  this  paper  we  first  describe  the 
BdG  formalism  for  local  defect  potentials  and  then  com¬ 
pare  with  the  new,  Koster-Slater  formalism.  Section  III 
describes  an  analytic  model,  bctsed  on  a  delta-function 
potential,  which  reproduces  some  of  the  quantitative  be¬ 
havior  of  the  numerical  results.  Section  IV  discusses  the 
results  of  the  numerical  calculations  for  magnetic  impu¬ 
rities,  non-magnetic  impurities,  impurities  incorporating 
both  magnetic  and  non-magnetic  potentials,  and  inhomo¬ 
geneities  in  the  pairing  interaction.  A  heuristic  picture  of 
the  electronic  structure  near  these  impurities  will  be  pre¬ 
sented  here,  and  the  calculations  will  be  compared  with 
the  analytic  model  of  Section  III.  Section  V  generalizes 
the  formalism  of  Section  II  to  the  case  of  strong-coupling 
and  anisotropically-paired  superconductors. 


II.  FORMALISM 


A.  Bogoliubov-de  Gennes  Equations 


To  place  our  new  formalism  in  context  we  will  contrast 
it  with  the  BdG  equations,  which  are  Schrodinger-like 
equations  for  the  electron  and  hole  components  of  the 
quasiparticle  wavefunctions  u(x)  and  v{x)  respectively. 
These  are,  for  a  free  electron  band  structure  with  mass 
m,  the  positive-energy  (E)  solutions  to 


2m 

2m 


-  E  +  Vo{x.)  +  (tVsCx) 

-  £  -  FoW  +  ctV5(x) 


Uo.(x)  +  A(x)i;o^(x)=  0 

+  A(x)ua(x)=  0. 

(1) 


Here  aVs  is  a  position-dependent,  spin-dependent  po¬ 
tential,  such  as  one  originating  from  an  impurity  with  a 
classical  spin.  Vo  is  a  position-dependent  non-magnetic 
potential  and  A(x)  is  the  inhomogeneous  order  parame¬ 
ter.  A(x)  can  be  chosen  real  since  the  defect  potential  is 
real®®.  The  quantization  direction  of  the  electronic  spins 
in  the  superconductor  {a  =  ±1/2)  is  chosen  parallel  to 
the  classical  spin.  A  classical  spin  has  no  quantum  dy¬ 
namics,  and  cannot  flip  the  quasiparticle  spin.  Hence 
spin  is  a  good  quantum  number  for  the  quasiparticles 
and  only  two  coupled  equations  (Eqs.  (1))  are  required. 
The  combinations  crVs  ±  Vq  have  physical  significance; 
aVs  +  Vo  is  the  potential  felt  by  an  electron  of  spin  a, 
while  aVs  -  Vq  is  the  potential  felt  by  a  hole  of  spin  a . 

The  spatially-dependent  order  parameter  is  deter¬ 
mined  self-consistently: 


2 


(7) 


A(x)  ^7(x,£:„„.)u„a(xX^(x)tanh  (2) 

na 

where  n  labels  the  states  for  each  spin  a,  T  is  the  temper¬ 
ature,  and  ks  is  Boltzmann’s  constant.  j{x,Encr)  is  the 
effective  electron-electron  interaction  potential,  which  is 

7(X,  Encr)~  lo  Ena  < 

=  0  Ena  >  (3) 

For  a  spherically-symmetric  defect  the  wavefunctions 
are  eigenstates  of  the  angular  momentum  with  quantum 
numbers  £  and  m.  Typically  the  defect  is  placed  in  a 
sphere  of  radius  R  with  appropriate  boundary  conditions. 
The  value  of  R  is  determined  by  the  spectral  resolution 
necessary  for  accurately  evaluating  Eq.  (2)  and  the  spec¬ 
tral  width  of  features  measurable  by  (for  example)  the 
STM.  The  typical  complications  resulting  from  approx¬ 
imating  an  infinite  system  by  a  finite-size  system  apply, 
such  as  discrete  states  above  the  energy  gap  and  the 
heavy  investment  of  computer  time  required  for  large 
values  of  i?.  For  example,  in  the  calculations  for  the 
vortex  in  NbSe2^'^’^^,  cf/^o  -  32  was  the  largest  ratio 
of  the  Fermi  energy  to  the  homogeneous  order  parame¬ 
ter  considered.  This  value  is  unrealistic,  and  is  a  result 
of  inappropriately  fitting  the  coherence  length  and  up¬ 
per  critical  field  of  NbSe2  to  a  free-electron  model.  A 
more  realistic  band  structure^®  has  a  bandwidth  to  en¬ 
ergy  gap  ratio  at  least  an  order  of  magnitude  greater, 
and  a  highly  anisotropic  Fermi  velocity.  A  system  where 
the  free-electron  model  is  appropriate,  such  as  niobium, 
has  an  cf/Aq  —  705  and  is  numerically  inaccessible  for 
general  potentials  Vs  and  Vo- 


B.  Koster- Slater  Inversion  Formalism 

We  now  introduce  a  self-consistent  method  which 
works  within  a  sphere  whose  radius  is  determined  by  the 
range  of  the  defect’s  potential  and  utilizes  the  contin¬ 
uum  spectrum  of  the  homogeneous  superconductor.  In 
essence,  we  invert  the  Gor’kov  equation  in  real  space. 
The  Gor’kov  equation"^^  for  a  defect  in  a  superconductor 
can  be  written  in  the  Nambu  formalism®^  as: 


F(x,x';cc;)- 

/CO 

'OO 

F*{x,x';u>)  = 

/OO 

dte^^9{t) {0|  {V’j  (x';  t) ,  (x;  0) }  jO) ,  (8) 

-OO 

-G;(x,x';-w)  = 

/OO 

dfe*‘^^^(t)(0|{i/>|(x';  t),^pl{x■,  0)}|0>.  (9) 

-OO 

The  explicit  subscripts  t  and  I  do  not  refer  to  the  spin 
of  the  excitation  in  the  superconductor  but  rather  to  the 
spin  band  of  the  normal-  state  used  to  construct  the  exci¬ 
tation .  The  key  concept  is  that  the  spin- -up  band  contains 
both  up  electrons  and  down  holes  just  as  the  spin-down 
band  contains  both  down  electrons  and  up  holes.  The 
convention  employed  here  is  standard  in  the  theory  of 
semiconductors  where  a  spin-up  electron  excited  above 
the  Fermi  energy  leaves  a  spin-down  hole  below  the  Fermi 
energy.  This  is  convenient  for  magnetic  potentials  since 
if  spin-up  electrons  are  attracted  to  a  magnetic  impurity 
spin-down  holes  should  be  repelled  by  the  impurity.  In 
the  presence  of  a  single  classical  impurity  spin,  the  quasi¬ 
particle  spin  is  a  good  quantum  number  despite  electron- 
hole  mixing.  Our  convention  determines  the  composition 
of  a  spin-up  quasiparticle  to  be  part  spin-up  electron  and 
part  spin-up  hole,  rather  than  part  spin-up  electron  and 
part  spin-down  hole. 

For  cj  >  0,  G(x,x';cj)  describes  spin-up  excitations, 
involving  the  mixing  of  electrons  in  the  spin-up  band 
with  holes  in  the  spin-down  band.  For  a;  <  0,  G(x,x';a;) 
describes  spin-down  excitations,  involving  the  mixing  of 
electrons  in  the  spin-down  band  with  holes  in  the  spin-up 
band.  Since  spin  is  a  good  quantum  number  it  is  not  nec¬ 
essary  to  use  a  4  X  4  formalism,  such  as  that  of  Ref.  62, 
The  notation  here  has  been  simplified  and  improved  rel¬ 
ative  to  Ref.  37. 

The  homogeneous  Green’s  function  g  is  independent 
of  cr,  so 


g(x,x'; 


g{x,x';ij) 

f*{x,x';u)) 


/(x,x';ai)  N 
-g*ix,x';-u;)  j  ' 


(10) 


j  dx"  [Six  -  x")  -  g(x,x";  w)V(x")]  G(x",x;  w) 

=  g(x,x';c4;)  (4) 

where  the  inhomogeneous  retarded  Green’s  function, 

/'‘^t(x,x';w)  Fix,x'-,u)  \ 
G(x,x,a;)  ^  f’*(x,x';a))  -G|(x,x';-w)  J'  ' 

The  elements  of  this  matrix  are 

Gt(x,x';w)=  • 

/OO 

dte^‘^^9(t)  (0|  {V't  (x';  t) ,  V-|  (x;  0)}  jO) ,  (6) 

-OO 


The  potential 


F5(x")  +  Vb(x") 
<SA(x") 


<5A(x")  \ 

Vsix")  -  Vo(x")  J 


(11) 


where  5A(x)  =  A(x)  -  Ao  and  V^s,  l^o,  and  A(x)  have 
similar  meaning  above  as  in  the  BdG  equations.  The 
factor  of  one-half  from  the  electron  spin  has  been  incor¬ 
porated  into  the  potential  Vs-  To  be  concrete,  and  with¬ 
out  loss  of  generality,  the  spin  direction  attracted  by  the 
potential  will  be  called  spin  up,  and  the  spin  direction 
repelled  will  be  called  spin  down®^.  The  self-consistency 
equation  for  the  order  parameter  is 
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/oo 

du)'y{'K,  a;)n(o;)ImF(x,  x;  cj)  -  (12) 

-oo 

where  n{uj)  is  the  Fermi  occupation  function  at  temper¬ 
ature  T, 

Self-consistent  spin-dependent  and  charge-dependent 
potentials  can  also  be  constructed  for  V(x)  using  the 
calculated  spatial  structure  of  the  spin  s(x)  and  charge 
p(x)  around  the  defect, 


s(x) 


Ini[Gt(x,x;cj)  -  Gj.(x,x;a))]\ 

TT  ) 

(13) 


ImGf7(x,  x;u;) 


(14) 


No  calculations  are  reported  for  such  potentials®^  in  this 
paper. 


For  a  BCS  superconductor  with  an  isotropic  gap  in  a  parabolic  band,  for  w  much  smaller  than  the  Fermi  energy, 


5(x,x';w) 

/(x,x';w) 


ttNq 

kpr 


^  y/Al  I  cos /c^rr  H- 


:  sin  kpr 


ttAoNo 


kpr^J  Al-uP' 


e  sinfc/rr 


(15) 


where  r  =  |x  —  x'|  and  No  is  the  density  of  states  for 
each  spin  at  the  Fermi  level.  The  coherence  length, 
^  =  TiUi^/TrAo,  where  vp  is  the  Fermi  velocity.  These 
expressions  are  valid  for  oj  above  and  below  A^  so  long 
as  the  imaginary  parts  of  both  /  and  g  are  multiplied  by 
sgna;. 

One  strength  of  this  formalism  is  its  reliance  on  the 
short-range  nature  of  the  defect’s  potential.  Solution  of 
Eq.  (4)  requires  inverting  the  frequency-dependent  real- 
space  matrix 

M(x,x';6.;)  =5(x-x')  -  g(x,x';  w)V(x').  (16) 


A  = 


Af  -Fn  A/^/ 


(17) 


where  n  and  /  label  the  near-field  region  and  far-field  re¬ 
gion,  respectively.  The  particular  example  of  M  is  block- 
triangular; 

I  -  0  \ 

_g/->nV  I ) 

M  g/-^nv(I  _  gn^nv)-l  I  )  * 


The  structure  of  M(x,x';a;)  allows  for  a  precise  descrip¬ 
tion  of  the  difference  between  the  near  field  and  the  far 
field.  We  require  that  the  defect’s  potential  V(x)  is  zero 
for  |x|  >  R.  The  space  |x|  >  R  belongs  to  the  far  field, 
whereas  the  space  |x|  <  R  belongs  to  the  near  field.  We 
can  then  separate  any  real-space  matrix  A  symbolically 
into  four  pieces: 


It  is  clear  from  Eq.  (18)  that  the  computational  effort  in 
inverting  M,  and  thus  finding  the  inhomogeneous  elec¬ 
tronic  structure,  is  entirely  determined  by  the  complexity 
of  inverting  We  also  obtain  the  useful  result  that 

the  electronic  structure  in  the  far  field  is  easily  deter¬ 
mined  once  the  electronic  structure  in  the  near  field  is 
known. 

The  local  density  of  states  in  the  near  field  is  merely 


-^^ImG,,(x,x;a;)  = -hm{23(x,x;w) -b  J  dx' [[(M”-"")  ®(x,x';w)]  g(x',x;w)] 

<T 

-jdx'  [(M"'^”)^“*^*(x,x';-a))  g*(x',x; | 
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while  in  the  far  field  it  is 


-  y^ImGg(x,x;cj)  =  -hm|25(x,x;  w) 

<T 

-b  j dx'dx"  [g(x,x';a;)V(x')  [(M"^")-'(x',x";w)]  g(x",x;a;)]„ 

-j  dx'dx"  [g*(x,x';-w)V(x')  [(M"^")(-')*(x',x";-u;)]  g*(x",x;  -a>)]^^  }. 


(19) 


(20) 
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In  the  limit  that  the  size  of  the  near-field  approaches  a  point,  Eq.  (20)  is 

--  ^ImG'a(x,x;w)  =  --Im|2p(x,x;w)  +  [g(x,0;w)Vefr(u;)g(0,x;  w)]jj  -  [g*(x,0;  -w)V*ff(-w)g*(0,x; 

a 

(21) 

The  inhomogeneous  part  of  the  right-hand  side  of  Eq.  (21)  has  the  following  form  if  the  off-diagonal  elements  of 
Vefr(w)  can  be  ignored  (this  is  appropriate  when  self-consistency  is  not  important): 

-ilm{([Veff(w)]ll  -  [Veff*(-w)]22)  (p(x,  0,  u;)ff  (0,  X,  Cc:)  -  /(x,  0,  w)/(0,  X,  w))}  .  (22) 

TT 


Hence  the  spatially-dependent  features  reported  in 
Ref.  26  are  qualitatively  retained  in  our  more  complete 
model.  For  a  magnetic  potential  in  the  Born  approxima¬ 
tion  there  would  be  no  such  signal,  since 


of  100  points  provides  a  numerically  robust  solution.  The 
self-consistent  solution  depends  on  the  value  of 

No^o/kp-  =  {2n^^kF)~^  =  (47r^€F/Ao)'"\  (26) 


Veff(a-)=(^^  (23) 

SO  the  distortions  in  the  LDOS  of  the  two  spin  directions 
are  equal  and  opposite. 

The  local  density  of  states  under  ideal  conditions  is  di¬ 
rectly  measured  by  a  scanning  tunneling  microscope.  The 
differential  conductivity  measured  at  a  point  x,  voltage 
V  and  temperature  T  can  be  related  to  the  local  density 
of  states  at  the  tip  location  as  follows: 

d/(x,y,T)  e2  ,  dnsTMiuj)^ 

(24) 

Here  e  is  the  charge  of  the  electron  and  nsTM{(^)  is  the 
Fermi  function  in  the  STM  tip 

(25) 

The  local  density  of  states  is  proportional  to  the  imagi¬ 
nary  part  of  the  retarded  Green’s  function  fully  dressed 
by  the  interaction  of  the  electronic  system  with  the  im¬ 
purity. 

Certain  features  of  the  above  equations  simplify  their 
numerical  implementation.  Each  angular  momentum 
channel  constitutes  an  independent  block-diagonal  sub¬ 
matrix  in  M(cj).  Since  the  bare  Green’s  functions  in 
Eqs.  (15)  have  analytic  expansions  in  spherical  harmon¬ 
ics,  M(cj)  can  be  calculated  quickly.  These  expansions 
are  detailed  in  the  Appendix. 

The  value  of  R  is  governed  by  the  longest-range  poten¬ 
tial.  In  this  paper  that  will  be  ^A(x).  We  model  Vs  and 
Vo  with  Gaussians  of  radius  and  evaluate  Eq.  (12)  at 
T  =  0.  The  (5A(x)  for  =  itNo\  f  dx\/5(x)|  =  1.75  and 
Vq  =  ttNoI/  dxyo(x)|  =  O.is  shown  in  Figure  1.  While 
oscillating  with  wavelength  ~  Trkp^,  5A(x)  falls  off  to  a 
negligible  potential  within  lOkp^.  A  typical  radial  grid 


nsTM((^)  = 


1  -h  exp 


/cj-eV 

I  kBT 


which  for  niobium  is  3.6  x  10”^  (Ref.  64).  This  is  the 
single  dimensionless  parameter  required  to  parametrize 
the  free-electron  model  of  a  superconductor. 

In  contrast  to  the  Koster-Slater  technique  described 
here,  which  exploits  a  physical  distinction  between  the 
near  field  and  the  far  field  to  accelerate  the  numerical 
calculation,  the  BdG  equations  treat  the  near-field  and 
the  far-field  on  the  same  level.  A  large  R  is  desired  for 
decent  spectral  resolution,  but  the  possible  size  of  R  is 
limited  by  numerical  constraints.  Hence  a  numerical  im¬ 
plementation  of  the  BdG  equations,  by  comparison,  is 
typically  slower  and  substantially  less  accurate  than  the 
Koster-Slater  technique. 


FIG.  1.  Change  in  the  local  order  parameter,  JA(x),  for 
a  magnetic  potential  strength  Vs  =  No  f  drVs(r)  =  1.75  as 
a  function  of  the  distance  from  the  impurity  r.  The  change 
becomes  negligible  beyond  lOkp^. 


III.  ANALYTIC  SOLUTION  OF  THE  POINT 
POTENTIAL 


A.  Magnetic  and  Non-Magnetic  Point  Potentials 
Approximating  the  local  potential  by  a  delta  function 
V(x)  =  V<5(x)  =  (  ^^  J  S(x)  (27) 

leads  to  a  simple  expression  for  Eq.  (4), 
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f  1  -5(0.0;w)(K5  +  Vo)  -f{0,0-uj){Vs  -  Vq) 
-f*{0,0;uj)iVs  +  Vo)  H-ff’(0,0;-w)(l/s-V'o) 


G(0, 0;  Lj)  -  M"-^"(w)G(0, 0;  w)  =  g(0, 0;  w). 


(28) 


In  principle  lVI(a;)  can  be  found  from  the  Green’s  func¬ 
tions  in  Eq.  (15),  however  there  is  a  divergence  in  the  real 
part  of  g{r]  cj)  as  r  ^  0.  This  divergence  is  coped  with  in 
Ref.  8  by  discarding  the  divergent  piece.  This  approxi¬ 
mation  is  essentially  an  assumption  of  strict  particle- hole 
symmetry  (not  merely  linearizing  c{k)  around  cf)-  We 
now  discuss  the  effects  of  this  approximation  on  the  lo¬ 
cal  properties  of  the  system  in  the  normal  state.  The 
heuristics  will  be  simpler  for  the  normal  state,  but  the 
conclusions  also  apply  to  the  superconducting  state. 

B.  Particle-Hole  Symmetry  in  the  Normal  State 

In  order  to  focus  on  the  spin-dependent  potential,  the 
non-magnetic  potential  will  be  set  to  zero.  The  nor¬ 
mal  state  properties  can  be  obtained  from  Eq.  (15)  for 
Ao  =  0.  That  yields  a  Green’s  function  appropriate  for 
an  outgoing  wave: 

=  (29) 

where  u)  is  considered  close  to  the  Fermi  surface  so  the 
change  in  momentum  due  to  a;  7^  0  is  negligible.  The 
inhomogeneous  local  density  of  states  for  a  delta-function 
potential  using  this  Green’s  function  is  unphysical,  since 
it  diverges  at  r  =  0.  Now  we  examine  the  local  density 
of  states  when  the  divergent  real  part  of  g  is  ignored. 

sin^  kFr\ 

(kpr)^  J 

(30) 

This  expression  yields  the  pathological  behavior  that  the 
local  density  of  states  near  a  spin-dependent  potential  is 
exactly  the  same  for  spin-up  electrons  as  for  spin-down 
electrons.  A  spatial  response  in  the  LDOS  to  a  spin- 
dependent  potential  that  is  identical  for  up  and  down  spin 
electrons  only  occurs  for  certain  band  structures  with  cf 
at  special  energies  (e.g.  half-filled  tight-binding  models). 
The  local  density  of  states  for  Vg  —  '^^o\ys\  ”  0.1  is 
plotted  in  Figure  2  as  a  function  of  distance  from  the 
potential.  Although  it  is  somewhat  distorted  from  its 
homogeneous  value,  it  does  not  show  the  spin-  depen¬ 
dent  asymmetry  of  a  more  realistic  potential.  An  exact 
calculation  for  a  Gaussian  of  range  kp^  is  also  shown  in 
Figure  2  for  comparison. 

A  more  realistic  approach  to  coping  with  the  diver¬ 
gence  in  Eq.  (29)  without  yielding  the  pathological  re¬ 
sult  of  Eq.  (30)  is  to  average  the  real  part  of  g  over  a 
range  given  by  the  eissiftned  range  of  the  potential.  That 
yields  a  finite  value  for  the  local  density  of  states  at  the 
potential,  but  does  not  control  the  behavior  for  small  r. 


To  perform  that  task,  we  consider  a  “muffin-tin”  Green’s 
function.  This  function  has  the  form 

r>Ro 

^  kFV 

=  —inNo  —  7ra  r  <  Rq  (31) 

where  a  is  the  average  of  the  real  part  of  g  over  the  range 
of  the  potential,  and  Ro  is  chosen  so  that  the  spatially- 
integrated  spectral  weight  of  the  Green’s  function  is  un¬ 
changed  (as  required  by  probability  conservation).  We 
show  in  Figure  2  the  local  density  of  states  calculated 
with  this  Green’s  function.  In  particular,  the  asymmetry 
in  the  response  of  spin-up  electrons  and  spin-down  elec¬ 
trons  to  a  spin-dependent  delta-function  potential  is  gov¬ 
erned  by  this  phenomenological  parameter  a.  The  agree¬ 
ment  with  the  exact  solution  is  good  at  the  origin  and 
far  from  the  impurity.  The  muffin-tin  Green’s  function  is 
discontinuous,  unfortunately,  but  yields  a  better  approx¬ 
imation  of  the  response  of  the  system  than  the  particle- 
hole  symmetric  approximation  (discarding  the  real  part 
of  the  Green’s  function). 

It  is  also  possible  to  generate  an  asymmetry  between 
spin-up  and  spin-down  electrons  by  adding  a  non- 


FIG.  2.  Spectral  weight  in  the  normal  state  around  a  mag¬ 
netic  impurity  with  =  0.1  as  a  function  of  the  distance  from 
the  impurity  r.  The  solid  lines  are  exact  solutions  for  spin-up 
(attracted  by  the  impurity)  and  spin-down  (repelled  by  the 
impurity)  electrons  for  a  Gaussian  potential  of  range  kp  in 
a  metal  with  a  free-electron  dispersion  relation.  The  dashed 
line,  which  is  the  same  for  both  spin-up  and  spin-down  elec¬ 
trons,  is  the  calculated  spectral  weight  for  the  particle-hole 
symmetric  delta-function  potential  model.  The  dot-dashed 
line  is  the  result  for  a  delta-function  potential  calculated  us¬ 
ing  the  muffin-tin  Green’s  functions.  The  muffin-tin  Green’s 
functions  fix  the  pathological  result  of  the  symmetric  model 
that  the  spectral  weight  is  the  same  for  spin-up  and  spin-down 
electrons.  The  muffin- tin  parameter,  a  =  0.704,  is  determined 
by  the  range  of  the  Gaussian  potential  and  is  therefore  not  a 
free  parameter. 


-^ImG<r(x,x;w)  =  No 


+ 


(t^NoVs)^ 


1 1  -  {wNoVs)^ 
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magnetic  potential  with  the  Vs  to  parametrize  the  im¬ 
purity,  but  still  maintaining  a  particle-hole  symmetric 
band  structure.  Since  the  response  of  a  particle-hole  sym¬ 
metric  system  to  a  potential  does  not  depend  on  the  sign 
of  that  potential,  the  vq  is  required  to  distinguish  between 
electrons  and  holes.  An  electron  feels  a  potential 
whereas  a  hole  feels  a  potential  Vg  -  vq.  This  additional 
non-magnetic  potential  only  breaks  particle-hole  symme¬ 
try  locally  (within  the  range  of  the  potential),  whereas 
for  a  realistic  band  structure  the  particle-hole  symmetry 
is  broken  everywhere  in  the  solid. 

For  a  Gaussian  potential  with  range  o, 


2  /  A  1  /  {kFa? 

~  ^kpa  V  2 


(32) 


For  the  Gaussian  potentials  numerically  calculated  in  this 
paper,  a  —  so  a  =  0.704. 


^(0,0;  w)  =  — ttA^o  + 


UJ 


(35) 


The  anomalous  Green’s  function  is  given  by  Eq.  (15), 
since  it  does  not  have  a  divergence  problem  as  r  — )•  0, 


f{0,0]u)  -ttNo 


(36) 


This  analytic  model  only  has  localized  states  in  the  £  = 
0  angular  momentum  channel,  as  expected  for  a  delta- 
function  potential.  Those  energies  are 


UJq 


VsSq  ±  [(2t;5(?o)^  -  4{v^  +  7^)((5o  -  7^)] 
1)2+72 


Ao 


(37) 


where 


C.  Self-consistency  within  the  Analytic  Model 


7  =  [(1  +  -  Vq  -  ^o)  “  2quo  -  1]  /2.  (38) 


As  seen  from  Figure  1,  the  distortion  of  the  order  pa¬ 
rameter  is  short-ranged  around  an  impurity.  We  may 
then  consider  the  effect  of  the  order-parameter  distortion 
on  the  electronic  structure  to  be  parametrized  by  an  effec¬ 
tive  delta-function  potential  (5AJ(x)  similarly  motivated 
to  the  delta-function  potentials  for  the  magnetic  and  non¬ 
magnetic  potentials.  The  potential  V(x)  is  changed  in 
the  following  way: 

V(x)  =  V<5(x)  =  (  )  J(x).  (33) 

The  relative  effect  of  the  SA  compared  to  the  other 
two  potentials  is  likely  to  be  small  for  the  potentials 
considered  in  this  paper.  Typically  NoVsfkp  ~  1  or 
NoVolkp  ^  1,  and  for  niobium  NoAo/k^  =  3.6  x 
10”^.  Even  for  a  small  coherence  length  of  ^ 

NoAojk’p  =  1.6  X  10“^.  For  convenience  we  define 
Sq  =  ttNoSA. 


Eq.  (37)  reduces  to  a  result  obtained  by  Shiba®  when 
=  a  =  ^0  =  0,  a  result  obtained  by  Rusinov®  when 
a  “  (5o  =  0,  a  result  obtained  by  Salkola,  Balatsky  and 
Schrieffer^®  when  a  =  0,  and  a  result  obtained  by  us^"^ 
when  uo  “  ^0  —  0* 

One  solution  of  Eq.  (37)  is  a  spin-up  quasiparticle  and 
the  other  is  a  spin-down  quasiparticle.  There  may  be 
only  one  real  solution  to  Eq,  (37);  then  only  one  £  =  0 
localized  state  exists  around  the  impurity.  This  occurs 
for  large  Vg,  When  =  0  the  localized  states  are  due 
to  order  parameter  suppression,  and  the  energies  of  the 
two  spin  states  are  degenerate.  This  follows  from  time- 
reversal  symmetry  in  the  absence  of  a  magnetic  potential. 
For  small  Vg  the  two  energies  are  split  by  an  amount 


Auo 


W2  +  72  • 


(39) 


D.  Energies  of  localized  states  in  the  Superconductor 


E.  Spectral  Weight  Asymmetry  in  the  Analytic 
Model 


The  energies  of  the  localized  states  of  angular  momen¬ 
tum  £  correspond  to  the  positive  energies  =  |n|,  where 

detM^^’^(n)  -  0,  (34) 

and  the  solution  is  traced  to  the  ^-channel  block  of  M 
(see  Appendix).  For  the  analytic  model,  M”^”(a;)  is 
the  matrix  shown  in  £q.  (28),  where  p(0,0;c(j)  in  the 
superconducting  statj^  is  tonstructed  similarly  to  that  of 
the  normal  state, 


A  spin-up  quasiparticle  consists  of  amplitudes  for  a 
spin-up  electron  (electron  in  a  spin-up  state) ,  and  a  spin- 
up  hole  (electron  missing  from  a  spin-down  state).  There¬ 
fore  the  spectral  weight  of  a  spin-up  localized  state  will 
be  divided  between  an  electron-like  pole  in  the  spin-up 
band  at  a;  =  (with  weight  A|(r;c4;o))  and  a  hole-like 
pole  in  the  spin-down  band  at  a;  =  — (with  weight 
A4,(r; -cjo)).  These  two  types  of  excitation  are  inde¬ 
pendently  resolvable  by  a  scanning  tunneling  microscope 
since  at  positive  sample  voltage  relative  to  the  tip,  the 
STM  places  electrons  in  the  sample,  whereas  at  negative 
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sample  voltage  the  STM  places  holes  in  the  sample.  We  define  the  energy  of  the  pole  in  the  spin-up  band  to  be  wt  and 
in  the  spin-down  band  to  be  w;.  Even  though  ojo  is  always  positive,  can  be  positive  or  negative,  and  uj^  =  -uj^. 

The  spatial  structure  of  the  spectral  weights  of  the  spin-up  band  and  spin-down  band  components  of  the  localized 
state  are  given  by 


,  /  ,  ttTVoAo  ,  \/A(,  oj 

Aa-  (r;  w)  =  — 5{(V  -  UJa) - 


2vs 


2  _  ,.,2  r -I-  (uo  -1-  Vs)uj‘^  +  (no  -I-  Vs)a^iAl  - 

A2 


-t- 


2(no  -t-  Vs)a  -  (1  -h  a'^){vl  -  v^)  (^q:(A2  -  +  w^A^  -i 


r  <  Ro 


3  2.  /  \Af^ 

i:NoAo  (Ag-a>^)i 
{kprY  2vsAl  ^ 

'  A2+a;2 


5{bj  —  Ua)  X 


’  A?  — 


-  ^0  +  (^5  - 'y  J  <  a - 


OJ 


sin^(A:F^) 


+fo)  1  +  (^s  " 


-  a  >  I  cos^(/cFr) 


U)^ 


-V2{Vs  +  Vo)  Vo  -  Vs  +  [1  -  q:(vo  -  Vs)] 


UJ 


UJ^ 


sm{kFr)  cos{kFr) 
r  >  Ro 


(40) 


The  above  expression  is  set  up  for  use  as  a  muffin-tin 
Green’s  function.  The  construction  of  such  a  Green’s 
function  requires  that  the  spectral  weight  integrates  to 
one.  However,  if  Ro  ->  0,  the  above  expression  can  be 
integrated  over  all  space,  and  still  yields  one-half  for  the 
spin-up  band  and  one-half  for  the  spin-down  band  to  or¬ 
der  Ao/ci?.  Hence  the  quasiparticle  is  one-half  electron 
and  one-half  hole.  This  provides  additional  confidence  in 
the  above  expressions. 

The  frequency-integrated  weight  at  the  defect  of  the 
two  types  of  excitation  can  be  calculated  within  the  an¬ 
alytic  model,  and  is  (for  (So  =  0) 

i4t(r  =  0)  _  1  -f  2a(uo  -  Vg)  +  (1  +  a^)(uo  -  Vs)^ 

^4,(r  =  0)  ”  1  -f  2a(uo  +  u*)  +  (1  +  O!^)ivo  +  Vs)^ 

For  a  spin-down  quasiparticle  there  is  an  electron-like 
pole  in  the  spin-down  band  and  a  hole-like  pole  in  the 
spin-up  band  and  the  relative  weight  is-still  given  by  the 
expression  above.  This  expression  for  vq  —  0  was  re¬ 
ported  in  Ref.  37. 

In  the  following  sec^tiori  these  results  are  compared  with 
the  numerical  calculations  of  properties  in  the  supercon¬ 
ducting  state  and  in  the  normal  state. 


r 


IV.  COMPARISON  WITH  NUMERICAL 
RESULTS  AND  DISCUSSION 

A.  Non-Magnet ic  Impurity 

Even  strong  non-magnetic  impurities  at  moderate  con¬ 
centrations  will  not  suppress  the  critical  temperature  of 
a  superconductor^.  Nevertheless,  it  was  recognized  early 
on^^  that  the  local  order  parameter  may  be  affected.  In 
Ref.  25  the  effect  of  the  non-magnetic  impurity  was  calcu¬ 
lated  in  the  far  field  by  modeling  the  impurity  potential 
with  a  phase  shift.  The  phase  shifts  were  evaluated  for 
two  models:  a  spherical  square- well  potential  and  a  delta- 
shell  potential.  Self-consistency  was  ignored  by  only  fo¬ 
cusing  on  regions  far  from  the  impurity  where  the  change 
in  the  order  parameter  is  small  compared  to  its  homoge¬ 
neous  value.  The  order  parameter  change  due  to  the  im¬ 
purity  was  found  to  oscillate  with  the  Fermi  wavelength, 
and  decay  as  r~^. 

Figure  3  shows  the  spectral  weight  A(r;cv)  at  several 
frequencies  above  the  energy  gap  near  a  strong  non¬ 
magnetic  impurity  with  a  Gaussian  potential  of  range 
k~^  and  strength  Vq  =  7/8,  calculated  self-consistently 
for  a  superconductor  with  ^k^  =  449  (niobium).  The 
spectral  weights  are  suppressed  to  approximately  30%  of 
their  homogeneous  value  at  the  center  of  the  potential. 
Only  continuum  states  are  shown  since  no  localized  states 
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FIG.  3.  Comparison  of  the  spectral  weight  in  the  super¬ 
conducting  state  at  several  frequencies  to  the  spectral  weight 
of  the  normal  state  and  the  spatial  dependence  of  the  or¬ 
der  parameter.  All  of  these  curves  are  indistinguishable,  but 
they  have  been  displaced  for  better  visibility.  The  spectral 
weights  are  normalized  to  their  value  at  large  distances  from 
the  impurity  r,  and  the  order  parameter  is  normalized  to  the 
homogeneous  order  parameter  Ao-  This  figure  indicates  that 
normal-state  properties  drive  the  spatially-dependent  features 
of  the  superconductor’s  spectrum. 


were  found  for  this  potential.  The  curves  showing  the 
spectral  weight  have  been  displaced  from  each  other  so 
that  they  may  be  distinguished.  Also  shown  displaced  in 
Fig.  3  is  the  spectral  weight  in  the  normal  state,  normal¬ 
ized  to  the  spectral  weight  in  the  homogeneous  metal. 
All  the  quantities  plotted  in  Fig.  3  are  identical  to  the 
accuracy  of  the  calculation.  Fig.  3  is  an  illustration  of 
a  relationship  between  the  spectral  weight  in  the  nor¬ 
mal  state  and  the  spectral  weight  in  the  superconducting 
state, 


A{r-,u;)  =  ^AM  (42) 

where  i4„(r)  is  the  spectral  weight  in  the  inhomoge¬ 
neous  normal  state  for  energies  near  the  Fermi  surface 
and  is  the  homogeneous  superconductor’s  spec¬ 

tral  weight  as  a  function  of  frequency.  2No  is  the  normal 
state’s  spectral  weight  far  from  the  impurity.  This  ex¬ 
pression  is  valid  for  small  r  and  small  cj,  the  regime  of 
interest  for  STM  on  a  superconductor.  For  uj  of  order 
Ao,  Eq.  (42)  is  valid  for  r  <  ^. 

We  further  illustrate  the  relationship  of  Eq.  (42)  in  Fig¬ 
ure  4(a),  which  shows  the  LDOS  for  this  non-magnetic 
impurity  as  a  function  of  voltage  and  position  calculated 
from  Eq.  (24)  with  T  =  —Ao/fcjg.  This  temperature 
corresponds  to  2K  for  niobium.  There  is  no  change  in 
the  energy  gap  due  to  this  non-magnetic  potential.  Fig¬ 
ure  4(a)  shows  that  it  is  merely  the  local  amplitude  of  the 
spectral  weight  which  is  reduced  —  this  would  manifest 
itself  in  a  locally  reduced  oscillator  strength  for  an  opti¬ 
cal  transition,  or  the  reduction  in  the  tunneling  current 
for  an  STM,  which  is  directly  proportional  to  the  LDOS 


shown  in  Figure  4(a).  Figure  4(b)  is  an  identical  calcu¬ 
lation  for  a  shorter  coherence  length,  ^  There 

appears  to  be  little  difference,  although  here  a  localized 
state  exists  near  the  continuum. 

Figures  3  and  4(ab)  show  for  the  non-magnetic  impu¬ 
rity  that  the  normal-state  electronic  structure  determines 
the  spatial  dependence  of  the  superconductor’s  A(r;a;) 
for  all  frequencies  including  near  the  energy  gap.  The 
potential  strength  of  the  impurity  is  orders  of  magnitude 
greater  than  Ao,  and  thus  locally  mixes  in  states  far 


FIG.  4.  Differential  conductance  (LDOS)  as  a  function  of 
voltage  and  position  calculated  around  a  non-magnetic  im¬ 
purity  with  Vo  =  0.875.  The  spectrum  is  suppressed  sub¬ 
stantially  in  the  vicinity  of  the  impurity.  The  temperature  is 
Ao/7.5kB,  which  for  niobium  corresponds  to  about  2K.  (a)  A 
coherence  length  appropriate  for  niobium,  ^kF  =  449.  (b)  A 
much  shorter  coherence  length,  ^kp  =  10. 
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from  the  Fermi  surface  in  the  homogeneous  metal.  These 
states  are  required  to  construct  probability  densities 
which  are  suppressed  by  70%  near  the  impurity.  The  spa¬ 
tial  structure  of  the  spectral  weight  in  the  normal  state 
is  essentially  identical  to  that  seen  in  Fig.  3  over  an  en¬ 
ergy  range  around  the  Fermi  surface  which  is  orders  of 
magnitude  greater  than  A^. 

Once  the  normal-state  band  structure  has  been  dis¬ 
torted  by  the  presence  of  the  non-magnetic  impurity, 
superconductivity  is  a  small  perturbation  within  a  nar¬ 
rower  frequency  range.  The  formation  of  the  gap  in  the 
single-particle  excitation  spectrum  in  the  superconduct¬ 
ing  state  is  characterized  by  the  “mixing”  of  electron 
and  hole  amplitudes  to  form  quasiparticles  near  the  gap 
edge.  These  quasiparticles,  therefore,  are  constructed 
from  single-particle  eigenstates  of  the  metal  which  have 
already  been  strongly  distorted  by  the  impurity  poten¬ 
tial. 

Equation  (42)  has  important  implications  for  spec¬ 
troscopy  on  a  superconductor,  for  one  of  the  procedures 
for  normalizing  spectra  taken  at  different  lateral  positions 
on  a  superconducting  surface  is  to  assume  that  the  LDOS 
at  a  particular  voltage  much  larger  than  Ao  is  the  same. 
This  is  an  attempt  to  correct  for  possible  changes  in  the 
tip-surface  distance  upon  moving  the  tip  laterally.  A 
small  change  in  the  tip-surface  distance  can  have  a  strong 
effect  on  the  tunneling  current.  Unfortunately  this  pro¬ 
cedure  will  prevent  an  experiment  from  seeing  changes  in 
the  LDOS  due  to  a  non-magnetic  impurity,  including  the 
conductance  oscillations  described  in  Ref.  26. 

We  now  discuss  the  properties  of  the  order  parame¬ 
ter.  A(r)  is  self-consistently  determined,  and  is  shown  in 
Fig.  3  for  small  r  to  be  identical  in  spatial  structure  to 
the  normal-state  spectral  weight, 

A(r)  ^  An(r)  ^  A{r-,uj) 

2No  Asciuj)'  ^  ’ 

Since  a  non-magnetic  potential  repulsive  to  electrons  at¬ 
tracts  holes,  and  A(r)  depends  equally  on  electron  and 
hole  amplitudes,  one  might  expect  a  non-magnetic  po¬ 
tential  to  have  little  effect  on  the  spatial  dependence  of 
the  order  parameter.  However  the  allowable  maximum 
spectral  density  of  holes  depends  on  the  spectral  density 
of  the  electron  band  where  the  holes  reside,  so  if  most 
electrons  are  excluded  from  the  site,  holes  will  be  effec¬ 
tively  excluded  as  well.  To  emphasize  this  point  we  note 
that  the  scaled  anomalous  spectral  weight  ImF(r,  r;  u)  is 
identical  to  the  scaled  A(r;c<;)  for  the  frequencies  shown 
in  Fig.  3  (and  for  all  relevant  frequencies  for  the  self- 
consistency  equation  Eq.  (12)).  Since  ImF(r,r;u;)  is 
proportional  to  the  product  of  electron  and  hole  ampli¬ 
tudes,  and  A(r;  uj)  is  proportional  to  the  electron  ampli¬ 
tude  squared,  the  spatial  structure  of  the  electron  and 
hole  spectral  weights  must  be  similar.  They  are,  since 
the  normal-state  speistral  weight  An{T]u)  is  roughly  fre¬ 
quency  independent  around  the  Fermi  energy  over  an  en¬ 
ergy  range  much  greater  than  Aq. 


We  now  comment  on  the  lack  of  localized  states  near 
the  non-magnetic  impurity  for  ^kp  =  449  and  the  small 
binding  energy  for  the  quasiparticle  for  ^kp  =  10.  The 
suppression  of  the  order  parameter  near  the  impurity 
may  be  considered  to  form  an  attractive  off-diagonal  po¬ 
tential  which  may  bind  quasiparticles.  Localized  states 
created  by  order-parameter  suppression  would  be  dou¬ 
bly  degenerate,  due  to  the  two  possible  spin  states  (see 
Eq.(37)).  Since  the  quasiparticle  is  half  hole  and  half 
electron,  if  the  electron  part  is  attracted  and  the  hole 
part  is  repelled,  one  might  expect  the  effects  of  such  a 
non-magnetic  potential  on  the  quasiparticle  to  cancel. 
However  the  binding  energy  of  the  localized  state  is  an 
order  of  magnitude  smaller  {u)o  =  (1  -  2  x  10""^)  A^)  than 
that  found  in  Section  IV. D  for  a  suppressed  order  pa¬ 
rameter  via  pairing  suppression.  This  may  be  explained 
by  the  well-known  repulsive  effect  (quantum  reflection) 
of  a  strong  attractive  potential  on  a  quantum-mechanical 
particle.  We  find  also  that  in  the  case  of  the  magnetic 
impurity  that  a  large  enough  non-magnetic  potential  of 
either  sign  will  suppress  the  binding  of  a  quasiparticle  to 
the  impurity.  We  note  here  that  the  ratio  of  the  non¬ 
magnetic  potential  to  the  off-diagonal  potential  (JA(r)) 
is  much  larger  (vq/NoAq  10^)  than  the  ratios  of  the 
non-magnetic  potentials  to  the  magnetic  potentials  con¬ 
sidered  below. 


B,  Magnetic  Impurity 

Recently,  we  have  presented  calculations  of  the  LDOS 
(and  thus  the  differential  conductivity  in  an  STM  exper¬ 
iment)  in  the  vicinity  of  a  magnetic  impurity^’^.  These 
calculations  indicated  that  the  spatial  structure  of  the 
electron  amplitude  of  the  localized  state  differed  strongly 
from  the  hole  amplitude  of  the  localized  state.  A  further 
result  was  that  the  spectrum  should  recover  to  the  ho¬ 
mogeneous  spectrum  within  a  few  atomic  spacings.  Sim¬ 
ilarly  motivated  calculations  of  the  LDOS  due  to  the 
i  =  0  localized  state  have  been  presented  since  then^®, 
although  these  calculations  did  not  address  the  contin¬ 
uum  LDOS.  The  two  models  used  in  Ref.  38  were  (1) 
a  J-function  model  solved  using  particle-hole  symmet¬ 
ric  Green’s  functions,  but  not  self-consistently,  and  (2) 
self-consistent  calculations  for  a  two-dimensional  tight- 
binding  s-wave  superconductor  within  the  BdG  equa¬ 
tions.  The  first  method  can  only  model  the  normal-state 
properties  properly  for  a  particle-hole  symmetric  band 
structure,  such  as  at  the  Van  Hove  singularity  in  a  two- 
dimensional  tight-binding  band  structure.  The  second 
method  must  contend  with  numerical  finite-size  effects, 
which  make  it  difficult  to  calculate  the  continuum  states. 
A  result  obtained  from  the  first  method  which  is  only  true 
for  special  band  structures  is  that  the  spatial  structure  of 
the  electron  and  hole  components  of  the  qu2isiparticle  are 
the  same.  The  authors  of  Ref.  38  did  raise  the  possibility 
of  an  additional  non-magnetic  potential  as  a  source 
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Energies  of  Quasiparticle  Poles 
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FIG.  5.  Energies  of  the  quasipaxticle  poles  as  a  function 
of  the  magnetic  potential  strength  Vs  for  the  first  two  an¬ 
gular  momentum  states  around  the  impurity.  The  spin-up 
and  spin-down  labels  refer  to  the  band  that  the  excitation 
resides  in  —  an  excitation  with  negative  energy  (hole-like) 
in  a  spin-up  band  is  a  spin-down  hole.  The  solid  line  cor¬ 
responds  to  =  10,  the  long-dashed  line  to  ^kp  =  449, 
the  short  dashed  line  to  the  symmetric  model  of  Shiba,  and 
the  dot-dashed  line  to  the  result  calculated  with  muffin-tin 
Green’s  functions  and  a  =  0.704.  At  a  critical  value  of 
Vs  =  v*Qj  the  up  poles  cross  to  negative  energies  and  the 
down  poles  cross  to  positive  energies,  indicating  a  change 
in  the  character  of"  the' ground  state.  The  kink  evident  in 
the  solid  and  long-dashed  lines  is  real,  and  due  to  the  dis¬ 
continuous  (at  T  =  0)  change  in  A(x)  at  u*o-  The  unim¬ 
portance  of  self-consistency  can  be  gauged  by  the  small  dif¬ 
ference  between  the  short-coherence  length  result  and  the 
long-coherence  length  result. 


of  electron-hole  amplitude  asymmetry  in  the  spatial 
structure  of  the  localized  state.  We  found^^  and  will 
explore  below  that  there  is,  for  realistic  band  structures, 
electron-hole  asymmetry  without  a  non-magnetic  poten¬ 
tial. 

We  will  begin  with  a  discussion  of  the  energies  of  the 
localized  states  around  a  magnetic  impurity  and  the  spin 
character  of  those  states.  Solutions  to  Eq.  (34)  can  be 
evaluated  numerically.  Figure  5  shows  the  dependence 
of  the  energies  of  the  localized  state  poles  for  the  first 
two  angular  momentum  channels  on  the  strength  of  the 
magnetic  potential.  Results  for  a  short  coherence  length 
(^  =  lOA:^^)  are  shown  (solid  line)  as  well  as  results  for 
a  long  coherence  length  (^  =  449A:^^).  The  localized 
quasiparticle  state  for  small  Vg  is  the  spin  state  attracted 
to  the  classical  spin,  which  we  will  label  up  (t)^^-  As 
the  potential  strength  increases,  the  excitation  energy  of 
each  angular  momentum  state  decreases.  At  some  critical 
value  the  localized  state  becomes  a  spin-down  exci¬ 
tation,  the  energy  changes  abruptly,  and  then  increases 
with  increasing  This  behavior  can  be  extracted  from 
the  analytic  model  (^q.  (37))  as  well. 

Also  shown  in  Figure  5  are  the  analytic  results  for  the 
pole  energies  for  a  —  0  (Ref.  8)  and  a  —  0.704.  The 


muffin-tin  model  is  no  better  than  the  particle-hole  sym¬ 
metric  model  in  predicting  the  localized  state  energies. 
The  muffin-tin  model  will  prove  more  successful  at  pre¬ 
dicting  the  asymmetry  between  the  electron  and  hole  am¬ 
plitudes  of  the  localized  state. 

The  unimportance  of  self-consistency  for  determining 
the  energies  of  the  localized  states  can  be  gauged  by  the 
small  difference  between  the  short-coherence  length  re¬ 
sult  and  the  long-coherence  length  result.  The  most  im¬ 
portant  feature  it  determines  in  Fig.  5  is  the  size  of  the 
discontinuity  in  the  localized  state  energy  at  This 
discontinuity  is  due  to  a  discontinuous  change  in  the  or¬ 
der  parameter  at  this  magnetic  potential  strength,  a  re¬ 
sult  pointed  out  in  Ref,  37  which  will  be  discussed  more 
below. 

A  magnetic  impurity  affects  each  member  of  a  pair 
oppositely.  An  up-spin  electron  is  attracted  to  the  im¬ 
purity  and  a  down-spin  electron  is  repelled.  The  differ¬ 
ence  in  energy  between  the  two  time-reversed  electron 
states  leads  to  a  locally  suppressed  pairing.  Therefore, 
the  energy  needed  to  break  a  pair  and  create  a  localized 
quasiparticle  with  angular  momentum  £  in  the  vicinity  of 
the  magnetic  impurity,  when  Vg  <  Vgi*,  is  reduced  from 
2Ao  to  Ao+iVi.  One  member  of  the  broken  pair  is  a  de¬ 
localized  spin-down  quasiparticle  at  the  gap  edge  (with 
energy  Ao).  The  other  member  of  the  broken  pair  is  a  lo¬ 
calized  spin-up  quasiparticle  with  energy  ue.  For  Vg  > 
there  is  a  spin-up  quasiparticle  with  angular  momentum 
£  in  the  ground  state  of  the  superconductor^^.  An  es¬ 
sential  point  about  the  ground  state  of  a  superconductor 
containing  classical  magnetic  impurities  is  that  at  T  =  0 
when  Vg  <  Vgt*  (for  all  £)  the  ground  state  is  composed 
entirely  of  paired  electrons.  When  Vg  >  Vgt*  for  any  £ 
then  the  T  =  0  ground  state  contains  localized  quasi¬ 
particles  as  well  as  pairs  and  there  are  new  low-energy 
excitations  including  the  re-formation  of  a  pair  as  well 
as  the  excitation  of  a  localized  quasiparticle  to  a  higher 
energy  localized  or  continuum  state. 

Evident  from  Eq.  (5)  is  that  the  electron  pole  involves 
entirely  single-particle  states  within  the  spin-up  band  for 
Vg  <  but  involves  entirely  states  within  the  spin- 
down  band  for  Vg  >  v*^.  The  hole-like  pole  always  in¬ 
volves  single-particle  states  in  the  opposite  spin  band 
from  the  electron-like  pole.  The  source  of  the  quasiparti¬ 
cle  amplitude  for  the  various  poles  is  indicated  in  Fig.  5. 

Figures  6(a),  6(b),  and  6(c)  show  LDOS  results  for 
Vg  =  0.52,  0.875,  and  1.75  respectively.  They  show  the 
state  split  off  from  the  continuum,  with  a  larger  electron¬ 
like  amplitude  than  hole-like  amplitude  (6(a)),  and  then 
lower  in  energy  with  an  increased  electron/hole  asym¬ 
metry  (6(b)).  Finally  the  larger  peak  becomes  hole-like 
(6(c)).  The  asymmetry  between  the  electron-like  and 
hole-like  peaks  becomes  more  pronounced  as  Vg  increases. 
We  note  that  the  larger  peak  is  always  associated  with 
the  spin-up  band,  whereas  the  other  is  associated  with 
the  spin-down  band.  Despite  the  apparent  differences 
in  peak  size,  the  spatially  integrated  electron  spectral 
weight  of  the  quasiparticle  is  equal  to  the  spatially  inte- 
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FIG.  6.  Differential  conductance  (LDOS)  calculated  around  a  magnetic  impurity  with  (a)  Va  =  0.5,  =  449,  (b)  Vs  =  0.875, 

^kp  =  449,  (c)  Vs  =  1.75,  ^kp  =  449,  (d)  Vs  =  0.875,  ^kp  =  10.  All  are  calculated  with  ksT  =  Ao/7.5-  In  the  progression  from 
(a)  to  (c)  the  asymmetry  between  the  two  peaks  increases  and  the  higher  peak  moves  to  lower  energies,  eventually  becoming 
hole-like. 


grated  hole  spectral  weight.  The  localized  quasiparticle 
is  always  half  electron  and  half  hole  for  all  potentials 
examined  here.  For  Vs  <  Vgo  spin-up  band  amplitude 
is  electron-like  and  the  spin-down  band  amplitude  is  hole¬ 
like.  At  v*Q  (1.32  for  niobium),  the  spin-up  component 
becomes  hole-like  and  the  spin-down  component  becomes 
electron-like,  as  required  by  the  change  in  the  spin  of  the 
excitation. 

Figure  6(d)  shows  the  LDt)S  for  a  markedly  different 
coherence  length,  ^  and  Vg  ~  0.875.  It  is  evalu¬ 

ated  for  the  same  values" of  AolksT  as  Figs.  6(abc)  and 
looks  almost  identical  to  Fig.  6(b).  Since  the  localized 
state  is  broadened  by  temperature  through  Eq.  (24),  this 


is  a  manifestation  of  the  proportionality  of  the  spectral 
weight  to  NoAo  (Eq.  (40)).  Figure  7  shows  the  spectral 
weight  at  the  origin  for  the  ^  =  0  state  and  for  the  £  =  1 
state  at  its  first  maximum  for  Vg  =  0.875  as  a  function  of 
the  inverse  of  the  coherence  length,  which  is  proportional 
to  NoAo-  Figure  8  shows  the  spectral  weight  for  ^  = 
449/cp^  as  a  function  of  Vg  for  the  spin-up  and  spin-down 
poles  at  the  origin  for  the  £  -  0  state  and  at  the  first 
maximum  for  the  £  =  1  state.  It  is  clear  that  a  non¬ 
magnetic  potential  is  not  necessary  to  obtain  an  electron- 
hole  asymmetry. 

In  Figure  9  we  show  the  asymmetry  at  the  impurity  as 
a  function  of  Vg  for  two  values  of  ^  —  a  long  coherence 


12 


1/(^kp) 


FIG.  7.  Spectral  weight  at  the  impurity  (r  =  0)  for  the  first 
angular  momentum  channel,  ^  =  0,  and  at  the  first  maximum 
for  the  second  angular  momentum  channel,  £  =  1,  for  poles 
in  both  the  spin-up  and  spin-down  bands,  as  a  function  of 
the  inverse  coherence  length,  showing  a  linear  behavior.  The 
magnetic  potential  strength  is  Vs  =  0.875. 

length  appropriate  for  niobium,  and  a  short  coherence 
length,  l^om  Figure  7  it  should  be  evident  that  the 
asymmetry  is  not  sensitive  to  It  is,  however,  pre¬ 
dicted  extremely  well  by  the  normal-state  spin-up  and 
spin-down  band  spectral  weight  asymmetry  at  the  impu¬ 
rity  (also  shown  in  Figure  9).  We  can  therefore  conclude 
that  as  with  the  non-magnetic  impurity  (Figure  4),  the 
spatial  structure  of  the  spectral  weight  is  a  normal-state 
property.  We  further  show  in  Figure  10  for  Vg  ~  0.875 
the  £  =  0,  spin-up  band  and  spin-down  band  projections 
of  the  normal-state  spectral  weight  to  compare  with  the 
localized  state  spin-up  band  and  spin-down  band  spec¬ 
tral  weights  for  two  values  oi  The  normal-state  and 
long-coherence  length  calculation  are  practically  indis¬ 
tinguishable.  The  insets  show  r^i4(r),  which  removes  the 
rapid  power-law  decay  of  the  state.  The  localized  states 
for  all  angular  momenta  £  will  decay  as  the  power  law 
r“^.  For  the  short-coherence  length  calculation  the 


FIG.  8.  Spectral  weight  at  the  impurity  (r  =  0)  for  the 
£  =  0  channel  for  =  449  ^  a  function  of  magnetic  po¬ 
tential  strength  for  poles  in  both  the  spin-up  and  spin-down 
bands.  The  spectral  weight  in  the  spin-up  band  pole  of  the 
£  =  0  localized  state  saturates  at  large  Vs .  Also  shown  are  the 
spectral  weights  at  the  first  maximum  for  the  £  =  1  localized 
states. 


Magnetic  Potential  Strength  (Vg) 


FIG.  9.  Ratio  of  the  spectral  weight  in  the  spin-up  band 
and  in  the  spin-down  band  at  the  impurity  (r  =  0)  as  a  func¬ 
tion  of  magnetic  potential  strength.  This  is  plotted  for  the 
normal  state  £  =  0  projected  Green’s  functions  (solid  line) 
as  well  as  for  the  localized  states  for  niobium  {^kp  =  449, 
long  dashed  line),  ioi  ^kp  —  10  (dot-dashed  line),  and  for  the 
muffin-tin  model.  The  muffin- tin  model  is  only  successful  for 
Vs  <  0.5,  but  that  is  due  to  a  breakdown  in  describing  the 
normal  state.  The  normal-state  electronic  structure  is  a  good 
predictor  of  the  superconductor’s  electronic  structure  for  the 
entire  range  of  Us. 

effect  of  an  exponential  envelope  is  also  visible.  In  the 
analytic  result  the  exponential  envelope  should  have  a 
range  R  =  7r^/2yi  -  (cjo/^oF?  which  corrects  to  better 
than  1%  the  discrepancy  in  Figure  10.  The  power-law 


FIG.  10.  Spectral  weights  for  £  =  0  localized  state  in  the  up 
and  down  bands  for  Vs  =  0.875.  The  solid  line  is  the  normal 
state  £  =  0  projected  spectral  weight,  the  long-dashed  line 
is  the  localized  state  in  a  superconductor  with  ^kp  =  449, 
and  the  dot-dashed  line  is  for  ^kp  =  10.  The  inset  shows 
the  spectral  weight  multiplied  by  to  remove  the  algebraic 
decay.  The  normal-state  and  long-coherence  length  results 
are  practically  indistinguishable.  The  deviation  shown  in  the 
short-coherence  length  superconductor’s  spectral  weight  is  fit 
to  within  1%  by  the  exponential  decay  factor  described  in  the 
text.  Hence  the  spatial  structure  of  the  spectral  weight  of 
the  superconductor’s  localized  state  is  well-predicted  by  the 
normal  state  spectral  weight. 
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fall-off  and  exponential  envelope  can  be  seen  directly 
from  Eq.  (15)  and  Eq.  (40). 

We  can  summarize  these  comments  with  a  general 
equation,  similar  in  concept  to  that  for  the  non-magnetic 
impurity,  Eq.  (42).  That  is,  for  a  localized  quasiparticle 
state  of  spin  the  spectral  weight-of-a  localized  state 
with  angular  momentum  i  would  be 

A„{T-,ui)  =  BAna (r, ^)e" ( ^ ^ )  <5(w  -  aa'ui), 

(44) 

where  5  is  a  normalization  factor  so  that  the  spectral 
weight  of  the  state  integrates  to  one,  and  is 

the  angular-momentum  i  projection  of  the  normal-si^Xe 
spectral  weight  in  the  spin  cr  band.  We  note  that  for 
small  r  there  is  an  approximate  relationship  between  the 
superconducting  state’s  spectral  weight  and  the  normal 
state’s  spectral  weight  in  each  spin  band, 

E 

dijA„{r-,u)  =  Ana{r),  (45) 

where  Ao  E  ^  This,  in  connection  with  Eq,  (44), 
implies  a  dependence  on  the  normal-state  structure  of  the 
continuum  spectral  weight  around  the  magnetic  impurity. 

We  now  return  to  the  structure  of  A(x).  This  quan¬ 
tity,  which  is  not  directly  observable,  formed  the  focus  of 
several  investigations  of  the  local  structure  around  a 


FIG.  11.  Order  parameters  as  a  function  of  distance  from 
the  impurity  r  calculated  for  several  magnetic  potential 
strengths  for  (a)  449  and  (b)  =  10.  In  both 

cases  there  is  a  discontinuous  change  in  the  order  parameter 
when  Vs  passes  through  the  critical  strength  vIq. 


FIG.  12.  Order  parameters  at  the  impurity  (r  =  0)  as  a 
function  of  magnetic  potential  strength  Vs  for  two  values  of 
the  coherence  length.  The  discontinuity  in  the  order  param¬ 
eter  at  dJo  is  much  larger  for  the  short-coherence-length  su¬ 
perconductor. 

magnetic  impurity.  The  oscillation  of  the  order  parame¬ 
ter  around  a  magnetic  impurity  was  first  evaluated  with¬ 
out  self-consistency^^“^^.  A  self-consistent  calculation  of 
the  order  parameter  at  the  impurity  and  very  far  away 
for  weak  impurity  potentials  was  done  by  Schlottmann^^ . 

As  shown  in  Figure  1,  for  large  values  of  A(x  = 
0)  <  0.  Sign  changes  in  A,  as  seen  in  pair  tunneling, 
have  been  suggested  for  magnetic  impurities  in  the  bar¬ 
riers  of  Josephson  junctions®®"®^.  The  sign  change  in 
A(0)  occurs  (at  T  =  0)  precisely  at  Due  to  the  spin 
and  frequency  symmetries  of  Eqs.  (4-11),  the  anomalous 
spectral  weight  ImF(r,  r,  uj)  associated  with  the  spin-up 
pole  is  always  equal  and  opposite  to  the  anomalous  spec¬ 
tral  weight  associated  with  the  spin-down  pole.  As  the 
pole  in  the  spin-up  band  goes  from  electron-like  {uj  >  0) 
to  hole-like  {u  <  0)  and  the  pole  in  the  spin-down  band 
goes  from  hole-like  to  electron-like  the  contribution  to 
A(0)  changes  sign  abruptly  at  T  =  0.  The  A(r)  resulting 
from  several  values  of  Vs  and  two  values  of  the  coherence 
length  are  shown  in  Figure  ll(ab).  The  discontinuity 
at  v*Q  is  more  pronounced  for  shorter  coherence  lengths 
since  the  localized  state’s  spectral  weight  is  more  con¬ 
centrated  at  the  impurity  (Eq.  (44)).  A(0)  as  a  function 
of  Vs  is  shown  in  Figure  12  for  two  values  of  the  coher¬ 
ence  length.  At  T  >  0  the  transition  would  be  smoothed 
somewhat. 

The  behavior  of  A(0)  as  a  function  of  Vg  comes  from 
the  introduction  at  v*q  of  a  quasiparticle  into  the  ground 
state  of  the  system.  The  spin-up  quasiparticle  localized 
near  the  impurity  in  the  ground  state  suppresses  the  lo¬ 
cal  order  parameter.  For  time-reversal  invariant  poten¬ 
tials  one  cannot  make  A(r)  negative  by  inserting  a  single 
quasiparticle,  since  the  suppression  from  one  quasiparti¬ 
cle  is  cancelled  by  the  lack  of  suppression  from  its  un¬ 
excited  Kramers  doublet  partner.  For  a  spin-dependent 
potential,  however,  the  anomalous  spectral  weight  near 
the  impurity  may  be  almost  entirely  contributed  by  the 
single  low-energy  localized  state.  When  a  quasiparticle  is 
present  in  the  ground  state,  the  ground  state  has  spin  | 
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up65,38  ^  negative  A(r)^’^.  Exciting  the  low-energy 

state  for  Vs  >  v*q  removes  the  spin-up  quasiparticle,  and 
therefore  increases  A(0),  whereas  excitation  of  quasi¬ 
particles  typically  reduces  A(x)  (which  is  the  case  for 
Vs  <  ‘^so)-  Also,  exciting  the  low-energy  state  reduces 
the  induced  spin  of  the  superconductor  at  the  impurity. 


C.  Combined  Magnetic  and  Non-Magnetic 
Potentials 

We  now  discuss  the  addition  of  a  non-magnetic  poten¬ 
tial  to  the  magnetic  potential.  It  had  been  suggested^® 
that  introducing  a  vq  with  a  Vs  would  provide  electron- 
hole  asymmetry.  We  find  that  it  does  change  the  asym¬ 
metry,  which  we  show  in  Figure  13  for  a  particular  Vs^  but 
that  once  again  this  is  a  normal-state  property.  The  rela¬ 
tionship  between  the  normal-state  spectral  weights  and 
the  superconducting-state  spectral  weights  of  Eq.  (44) 
still  holds.  Introducing  vq  also  alters  the  localized-state 
energies  (see  Eq.  (37)),  which  we  show  in  Figure  14  for 
Vs  =  0.875  and  —  449.  The  presence  of  a  non¬ 
magnetic  potential  may  affect  the  value  of  We 

show  in  Figure  15  a  partial  diagram  of  the  ground  state  as 
a  function  of  the  parameters  Vs  and  vq  for  (kp  =  449.  We 
note  that  the  boundary  between  the  two  ground  states  is 
not  shifted  much  from  ^kp  =  449  to  ^kp  =  10,  hence  the 
condensation  energy  is  not  very  significant  in  determin¬ 
ing  this  boundary.  This  fact,  and  the  observation  that 
the  energy  of  the  localized  state  (in  the  absence  of  self- 
consistency)  relative  to  the  homogeneous  order  parame¬ 
ter  does  not  depend  on  any  parameter  of  the  supercon¬ 
ducting  state  (Eq.  (37))  suggests  that  an  understanding 
of  the  transition  at  v*o  vaecy  be  assisted  by  normal-state 
concepts.  We  hope  to  address  this  issue  in  a  future 


Non-Magnetic  Potential  Strength  (v^) 


FIG.  13.  Ratio  of  the  spectral  weight  at  the  impurity  in 
the  spin-up  band  to  the  spin-down  band  for  Vs  =  0.875  and 
=  449  as  the  non-magnetic  potential  vo  varies.  In  a 
similar  result  to  that  seen  in  Fig.  9,  the  normal-state  spectral 
weights  are  good  predictors  of  the  superconducting  state’s 
spectral  weight.  We  note  that  the  curves  are  not  symmetrical 
around  vo  =  0,  which  iesults  from  a  realistic  band  structure 
without  peirticle-hole  symmetry. 


Non-magnetic  Potential  Strength  (Vq) 


FIG.  14.  The  energy  of  the  spin-up  pole  is  shown  as  a  func¬ 
tion  of  non-magnetic  potential  strength  for  Vs  =  0.875  and 
^kp  =  449.  The  energy  of  the  spin-down  pole  is  just  the 
negative  of  the  energy  of  the  spin-up  pole.  The  analytic  mod¬ 
els  do  not  perform  particularly  well  in  reproducing  the  pole 
energies,  although  the  muffin-tin  model  does  show  a  simil2Lr 
asymmetry  around  uo  =  0  to  the  numerical  calculations. 


publication. 

We  now  mention  an  example  designed  to  showcase 
the  drawbacks  of  Ginzburg-Landau  theory.  We  present 
in  Figure  16  the  order  parameter  for  a  mixed  magnetic 
and  non-magnetic  impurity  (ug  =  7/8,  vq  =  —7/4,  and 
^kp  =  449).  The  order  parameter  is  everywhere  larger 
than  in  the  homogeneous  superconductor,  however,  the 
presence  of  a  localized  state  within  the  gap  indicates  that 
superconductivity  has  been  weakened  around  the  impu¬ 
rity.  Since  Ginzburg-Landau  theory  focuses  on  the  or¬ 
der  parameter  and  ignores  the  quasiparticle  structure,  a 
Ginzburg-Landau  perspective  would  incorrectly  predict 
an  enhancement  of  superconductivity  in  the  region. 

We  conclude  by  attempting  to  provide  some  guidance 
for  attempts  to  extract  impurity  potentials  from  STM 
measurements.  Figure  17  shows  normal-state  dl/dV^s  for 
various  potentials.  These  curves  should  also  represent 


Localized  Quasiparticle 
in  the  Ground  State 


No  Localized  Quasiparticle 
in  the  Ground  State 

-1  ’  0  '  ‘“'t 
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FIG.  15.  Calculated  boundary  between  two  ground  states 
around  the  magnetic  impurity  for  ^kp  —  449.  For  a  large 
enough  magnetic  impurity  strength  a  quasiparticle  is  bound  in 
the  ground  state.  The  minimum  magnetic  impurity  strength 
depends  on  the  non-magnetic  impurity  strength.  For  still 
larger  magnetic  impurity  strengths  there  would  be  ground 
states  with  more  than  one  bound  quasiparticle. 
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FIG.  16.  The  change  in  order  parameter  axound  an  attrac¬ 
tive  non-magnetic  potential  combined  with  a  magnetic  po¬ 
tential  for  =  449.  The  order  parameter  is  larger  at  the 
impurity  than  in  the  homogeneous  superconductor,  yet  there 
exists  a  localized  state  in  the  gap. 

the  frequency-averaged  spectral  weight  measured  in  the 
superconducting  state  (see  Eq.  (45))®^. 

The  enhancement  or  suppression  of  spectral  weight 
near  the  origin  is  particularly  sensitive  to  Uq  -  A  measure¬ 
ment  of  this  quantity,  the  energy  of  the  localized  state 
and  the  asymmetry  of  the  electron  and  hole  amplitudes  at 
the  impurity  overconstrains  vq  and  Vs ,  given  an  assump¬ 
tion  of  the  shape  of  the  potential.  To  extract  informa¬ 
tion  about  the  potential’s  detailed  shape  would  require 
a  fitting  procedure  using  the  differential  conductivity  at 
various  positions.  If,  for  some  reason,  the  spin-down  am¬ 
plitude  were  too  small  to  measure,  it  may  remain  possible 
to  constrain  the  potential  strength  using  the  frequency- 
integrated  spectral  weight  and  the  localized-state  energy. 

It  seems  appropriate  to  mention  again  the  tendency 
to  normalize  spectra  according  to  the  LDOS  measured 
at  energies  much  larger  than  Ao.  Since  the  normal- 
state  LDOS  near  the  impurity  changes  substantially  in 
the  presence  of  magnetic  or  non-magnetic  impurities,  an 
experiment  performed  using  such  a  normalization  pro¬ 
cedure  would  yield  impurity  parameters  of  questionable 
validity. 

D.  Pairing  Suppression 

The  pairing  potential,  7(x)  in  Eq.  (12),  may  also  have 
spatial  structure.  When  this  parameter  is  changed  it  in¬ 
duces  a  change  in  the  order  parameter  which  produces 
an  off-diagonal  potential  felt  by  the  qucisiparticles.  We 
set  Vs  Vq  =  0  so  that  there  is  no  magnetic  or  non¬ 
magnetic  potential  to  compete  with  the  order  parameter 
change.  Figure  18  shows  the  order  parameter  around  a 
short-range  suppression, 

7(x)  =  [l  -  7o,  (46) 

for  two  values  of  the  coherence  length.  The  order  param¬ 
eter  is  strongly  suppressed  and  since  7(0)  =  0,  A(0)  =  0. 
For  long  coherence  lengths,  this  change  in  the  order  pa¬ 
rameter  hss  no  effect  on  the  local  density  of  states  (shown 
in  Figure  19(a)).  Itds  possible  to  localize  quasiparticle 
states,  however,  at  shorter  coherence  lengths.  These 


r  (units  of  kp"^) 


FIG.  17.  Spectral  weights  in  the  normal  state  as  functions 
of  the  distance  from  the  impurity  for  several  combinations 
of  magnetic  and  non-magnetic  potentials  (a)  only  a  magnetic 
potential,  and  (b)  a  fixed  magnetic  potential  =  0.875  and 
a  varying  non-magnetic  potential.  By  making  measurements 
around  the  impurity  in  the  normal  state  (or  integrating  the 
superconductor’s  spectrum  over  a  frequency  much  larger  than 
Ao),  information  about  the  structure  of  the  impurity  may  be 
obtained. 

can  produce  features  in  the  local  density  of  states  which 
are  visible.  One  such  case  is  shown  in  Figure  19(b).  The 
energy  of  the  localized  state  is  cuo  =  (1  —  4  x  10~^)Ao. 
Whereas  a  non-magnetic  potential  changes  the  local 


r  (units  of  kp  ^) 

FIG.  18.  Order  parameter  as  a  function  of  distance  r  from  a 
defect  with  a  suppressed  pair  potential,  but  no  single-particle 
potential. 
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FIG.  19.  Differential  conductance  (LDOS)  around  a  defect  with  a  suppressed  pair  potential,  but  no  single-particle  potential, 
(a)  ^kp  =  449,  ksT  =  Ao/7.5  —  there  is  no  evidence  of  any  change  in  the  spectrum  due  to  the  order  parameter  suppression 
shown  in  Fig.  18.  (b)  ^kp  —  10,  ksT  —  Ao/lOO  —  a  localized  state  very  near  the  continuum  enhances  the  continuum  edge 
seen  in  tunneling  near  the  defect. 

density  of  states  without  significantly  changing  the  energy  gap,  a  pairing  suppression  has  a  very  weak  effect  on  both, 
especially  in  the  long-coherence  length  limit. 


V.  STRONG-COUPLING  AND  ANISOTROPIC  ORDER  PARAMETERS 


A  few  observations  are  in  order  concerning  the  extension  of  this  formalism  to  systems  where  the  homogeneous 
order  parameter  has  important  frequency  or  spatial  structure.  The  Gor’kov  equation  (Eq.  (4))  changes  due  to  the 
more  general  form  for  the  off-diagonal  potential  originating  from  the  order  parameter.  Taking  this  opportunity  to 
generalize. 


j  dx"dx"'  [^(x  -  x'')S{x"  -  x'")  -  g{x,x"-uj)Yix",x"'-,io)]  G(x"',x';w)  =  g(x,x';w) 


where 


_  /  -  x'")  <5A(x'',x"» 


V(x,x,w)-I  5A(x",x"';a;)  Vh^{x";cj)6{x"  -  x'") 


)■ 


(47) 


(48) 


The  diagonal  terms  are  general  potentials,  possibly 
frequency-dependent,  effective  on  spin-up  electrons  (14t) 
and  spin-up  holes  (Vh^)  for  a;  >  0.  Since  this  poten¬ 
tial  is  diagonal  in  frequency,  as  is  the  Gor’kov  equation, 
the  frequency  structure  of  the  order  parameter  does  not 
add  any  additional  complication  to  numerically  solving 
the  Gor’kov  equation.  However,  the  addition  of  spatial 
structure  to  the  pairing  has  added  another  integral  over 
the  volume  to  the  Gor’kov  equation,  and  thus  dramati¬ 
cally  increased  the  size  of  the  matrix  which  needs  to  be 
inverted.  Fortunately,  the  range  of  the  order  parameter 
in  |x  -  x'l  is  truncated  by  the  range  of  the  pairing  inter¬ 
action.  For  the  isotropic-gap  superconductor  we  have 


considered  for  most  of  this  paper,  the  effective  pairing  in¬ 
teraction  is  modeled  by  a  delta  function  in  space.  It  is 
possible  to  obtain  anisotropic  order  parameters,  includ¬ 
ing  the  d-wave  order  parameters  possibly  appropriate  for 
high-temperature  superconductors,  merely  by  consider¬ 
ing  pairing  with  nearest  neighbors  on  a  square  tight- 
binding  lattice.  This  model  has  been  implemented  with 
the  BdG  formalism  for  non-magnetic  impurities^^"^^  and 
magnetic  impurities^®  and  vortices^^  in  d-wave  supercon¬ 
ductors.  When  formulated  on  a  lattice  the  addition  of 
nearest-neighbor  pairing  multiplies  the  rank  of  the  matrix 
by  1  -h  z,  where  z  is  the  coordination  number 
of  the  lattice. 
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The  order  parameter’s  frequency  dependence  compli¬ 
cates  the  self-consistency  equation  (Eq.  (12)).  It  must 
now  be  solved  for  each  frequency: 


A(r,r';cj)Z(r,r';w) 

/(» 

den{e)  [ReF(r,  r';  e)]  {K+  (r,  r' ;  e,  w)  -  t/c) 

-OO 

w[l-Z(r,r';w)] 

/OO 

den(e)  [ReG(r,r'; e)]  R_(r,r';e,t<;),  (49) 

-OO 


where  K±  are  kernals  of  the  pairing  interaction  and  are 
different  for  each  mechanism  of  superconductivity.  They 
can  be  determined  from  the  homogeneous  solution.  Uc 
is  a  Coulomb  factor.  Z,  the  quasiparticle  weight,  is 
solved  for  self-consistently.  Incorporating  these  strong 
coupling  effects  allows  a  determination  of  the  effect  of 
the  frequency-dependence  of  the  pairing  interaction  on 
the  electronic  structure  around  a  defect.  While  selecting 
the  particular  location  where  the  STM  tip  is  similar  to 
selecting  the  momentum  of  the  quasiparticles  of  interest, 
selecting  the  STM  voltage  indicates  which  order  param¬ 
eter  frequency  one  wishes  to  probe. 


the  inhomogeneous  normal  metal  multiplied  by  a  decay¬ 
ing  exponential  whose  range  is  determined  by  the  energy 
of  the  localized  state.  Since  the  spin-up  band  LDOS  in 
the  normal  state  differs  from  the  spin-down  band  LDOS 
in  the  normal  state,  the  electron-like  pole  of  the  localized 
quasiparticle  will  have  different  spatial  structure  than  the 
hole-like  pole  of  the  localized  quasiparticle. 

The  self-consistent  calculations  described  here  have 
been  performed  with  a  new,  powerful,  Koster-Slater  tech¬ 
nique  which  allows  the  Gor’kov  equation  to  be  solved  in 
principle  exactly.  Although  we  have  only  presented  cal¬ 
culations  for  weak-coupling  isotropic  order  parameters 
within  a  free-electron  model,  we  have  formulated  the  ex¬ 
tension  of  the  new  technique  to  strong-coupling  pairing 
and  general  band  structures  and  order  parameter  sym¬ 
metries. 
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VI.  SUMMARY 

The  local  electronic  structure  around  an  defect  reflects 
both  the  properties  of  the  defect  and  the  medium  it  is 
embedded  in.  For  a  strong  non-magnetic  or  magnetic  im¬ 
purity  in  a  superconductor,  the  distortion  of  the  normal- 
state  properties  by  the  strong  impurity  plays  a  vital  role 
in  the  response  of  the  superconducting  medium.  The 
LDOS  for  the  inhomogeneous  superconductor  can  be  re¬ 
lated  to  the  LDOS  for  the  inhomogeneous  normal  metal 
via  equations  like  Eq.  (42)  and  Eq,  (44).  In  the  case  of 
a  non-magnetic  impurity  with  no  localized  states  around 
it,  the  LDOS  for  the  inhomogeneous  superconductor  is 
merely  the  normalized  LDOS  for  the  inhomogeneous  nor¬ 
mal  metal  multiplied  by  the  homogeneous  superconduc¬ 
tor’s  density  of  states.  This  should  suggest  some  caution 
for  the  method  used  to  normalize  STM  spectra  taken  at 
different  places  on  a  superconductor’s  surface. 

For  the  case  of  a  localized  state  (such  as  around  a  mag¬ 
netic  impurity)  with  angular  momentum  the  LDOS  for 
the  state  is  the  angular-moment um-^-projected  LDOS  of 


APPENDIX 

The  expansion  of  the  Green’s  functions  of  Eq.  (15) 
suitable  for  a  three-dimensional  spherically  symmetric 
situation  are  detailed  here.  The  homogeneous  Green’s 
functions  depend  on  r  and  r'  through  r,  r',  and  cos  7  = 
(r  •r')/rr'.  Then  the  Green’s  functions  can  be  expanded 
in  Legendre  polynominals  Pi{cos^). 

g(r,r';w)= 

9^  I  1 

g{r,  r',cos7;w)=  —^'^9e{r,r']uj)Pe{cos'Y)  (50) 

I 

and 

ge{r, r' ■, uj)  =  2Tr  j  Pi{x)g{r,r' ,x]u)).  (51) 

Evaluating  Eq.  (51)  for  both  the  normal  and  anomalous 
Green’s  functions  of  Eq.  (15)  yields 


5,(r,r»  =  [i  (l  +  ;^^) 

i 


ft{r,r'-,(J)  =  - 


V'w2  _  1 

3  1 


J,+i({l  -  ^^2  _  l/e}r<)/f2i({l  -  n/u;2  _ 


(52) 


J,+  i({l  +  -  l/^}r>) 
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(53) 


+Je+i{{l  -  n/u.2  -  -  l/Or>) 


where  Jt,  and  are  standard  Bessel  functions, 
(r>)  is  the  smaller  (larger)  of  r  and  a;  is  in  units 
of  Ao  and  r  and  r'  are  in  units  of  kp^.  The  Green’s 
functions  are  in  units  of  A^o- 

The  Gor’kov  equation,  Eq.  (4),  can  now  be  written  in 
a  form  diagonal  in 


Gi{r,r']uj)=  (54) 

rOO 

ge{r,r';u)+  /  r^cir„gf(r,r„;  w)V(r„)G^(r„,r';a;). 

Jo 


Thus  the  three-dimensional  integral  has  been  reduced  to 
a  one-dimensional  radial  integral.  Since  the  numerical  in¬ 
version  procedure  depends  on  inverting  a  matrix  whose 
rank  is  proportional  to  the  number  of  spatial  points  con¬ 
sidered,  this  reduction  to  a  one- dimensional  integral  dra¬ 
matically  increases  the  speed  of  this  calculation  over  a 
calculation  for  a  three-dimensional  potential  which  is  not 
spherically  symmetric. 
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